diff --git a/docs/user_guide/how_to_interface_with_user_defined_model.md b/docs/user_guide/how_to_interface_with_user_defined_model.md index 9f13197f5..173de4ed3 100644 --- a/docs/user_guide/how_to_interface_with_user_defined_model.md +++ b/docs/user_guide/how_to_interface_with_user_defined_model.md @@ -25,8 +25,8 @@ You can combine an existing H2Integrate model and a custom model for the same te To demonstrate this capability, we include a minimal example of a custom technology model: a **paper mill**. This example includes: -- A `PaperMillPerformance` model that converts electricity input to paper output. -- A `PaperMillCost` model that estimates capital and operational expenditures. +- A `CustomPaperMillPerformance` model that converts electricity input to paper output. +- A `CustomPaperMillCost` model that estimates capital and operational expenditures. - A `PaperMillFinance` technology finance model that computes the levelized cost of paper production (LCOP). Refer to the [Paper Mill Model Example](https://github.com/NatLabRockies/H2Integrate/tree/develop/examples/06_custom_tech/) for a complete walkthrough. diff --git a/examples/06_custom_tech/tech_config.yaml b/examples/06_custom_tech/tech_config.yaml index 62bb1b6ae..d2f200661 100644 --- a/examples/06_custom_tech/tech_config.yaml +++ b/examples/06_custom_tech/tech_config.yaml @@ -30,10 +30,10 @@ technologies: cost_year: 2019 paper_mill: performance_model: - model: PaperMillPerformance + model: CustomPaperMillPerformance model_location: user_defined_model/paper_mill.py cost_model: - model: PaperMillCost + model: CustomPaperMillCost model_location: user_defined_model/paper_mill.py finance_model: model: PaperMillFinance diff --git a/examples/06_custom_tech/user_defined_model/paper_mill.py b/examples/06_custom_tech/user_defined_model/paper_mill.py index f6ef47c7e..a353c1b8b 100644 --- a/examples/06_custom_tech/user_defined_model/paper_mill.py +++ b/examples/06_custom_tech/user_defined_model/paper_mill.py @@ -15,7 +15,7 @@ class PaperMillConfig(BaseConfig): electricity_usage_rate: float = field() -class PaperMillPerformance(om.ExplicitComponent): +class CustomPaperMillPerformance(om.ExplicitComponent): _time_step_bounds = ( 3600, 3600, @@ -57,7 +57,7 @@ class PaperMillCostConfig(CostModelBaseConfig): plant_capacity: float = field() -class PaperMillCost(CostModelBaseClass): +class CustomPaperMillCost(CostModelBaseClass): _time_step_bounds = ( 3600, 3600, diff --git a/examples/36_paper_mill/36_paper_mill_mn.yaml b/examples/36_paper_mill/36_paper_mill_mn.yaml new file mode 100644 index 000000000..484131214 --- /dev/null +++ b/examples/36_paper_mill/36_paper_mill_mn.yaml @@ -0,0 +1,6 @@ +name: H2Integrate_config +system_summary: This reference paper mill plant is located in Minnesota and for its first pass, it contains paper mill plant + powered by grid. The system is designed to produce paper at a constant rate throughout the year. +driver_config: driver_config.yaml +technology_config: tech_config.yaml +plant_config: plant_config.yaml diff --git a/examples/36_paper_mill/Breakdown_cost_plots.py b/examples/36_paper_mill/Breakdown_cost_plots.py new file mode 100644 index 000000000..26fc39412 --- /dev/null +++ b/examples/36_paper_mill/Breakdown_cost_plots.py @@ -0,0 +1,185 @@ +# """ +# Created on Fri May 15 07:38:06 2026 + +# @author: mkoleva +# """ + +# import pandas as pd +# import matplotlib.pyplot as plt + +# # Load Excel file +# file_path = "Breakdown_costs_per_scenario.xlsx" +# df = pd.read_excel(file_path, sheet_name="Sheet1", header=None) + +# # Scenario labels — edit however you prefer +# scenarios = [ +# "Paper + Pulp", +# "SAF with H2", +# "SAF with low-carbon H2", +# "Paper + Pulp\nSAF with H2", +# "Paper + Pulp\nSAF with low-carbon H2" +# ] + +# # Extract cost component names +# components = df.iloc[3:, 0].values + +# # Build scenario value arrays (sum of appropriate columns) +# records = {} +# records["Paper + Pulp"] = df.iloc[3:, [1, 2, 3]].astype(float).sum(axis=1).values +# records["SAF with H2"] = df.iloc[3:, [3]].astype(float).sum(axis=1).values +# records["SAF with low-carbon H2"] = df.iloc[3:, [4]].astype(float).sum(axis=1).values +# records["Paper + Pulp\nSAF with H2"] = df.iloc[3:, [5, 6, 7]].astype(float).sum(axis=1).values +# results = df.iloc[3:, [8, 9, 10]] +# records["Paper + Pulp\nSAF with low-carbon H2"] = results.astype(float).sum(axis=1).values + +# # Build DataFrame +# plot_df = pd.DataFrame(records, index=components) +# plot_df = plot_df[scenarios] # order consistently + +# # Assign custom colors +# colors = [] +# for comp in plot_df.index: +# if "CapEx" in comp: +# colors.append("navy") +# elif "OpEx" in comp: +# colors.append("orange") +# elif "Feedstock" in comp: +# colors.append("deepskyblue") +# elif "Taxes" in comp: +# colors.append("lightpink") +# elif "Finances" in comp: +# colors.append("yellowgreen") +# else: +# colors.append(None) # Let matplotlib choose default + +# # Plotting +# plt.figure(figsize=(10, 6)) +# bottom = [0] * len(scenarios) + +# for idx, comp in enumerate(plot_df.index): +# plt.bar( +# scenarios, +# plot_df.loc[comp], +# bottom=bottom, +# color=colors[idx], +# label=comp +# ) +# bottom = [bottom[i] + plot_df.loc[comp][i] for i in range(len(scenarios))] + +# plt.xlabel("Scenario") +# plt.ylabel("Cost ($/kg)") +# plt.title("Cost Breakdown per Scenario") + +# # FORCE horizontal x-axis labels +# plt.xticks(rotation=0, ha="center") + +# plt.legend() +# plt.tight_layout() + +# plt.savefig("stacked_cost_breakdown_final.png", dpi=300) +# plt.show() + +import textwrap + +import numpy as np +import pandas as pd +import matplotlib.pyplot as plt + + +# ------------------------------------------------------------- +# LOAD EXCEL +# ------------------------------------------------------------- +file_path = "Breakdown_costs_per_scenario.xlsx" +df = pd.read_excel(file_path, sheet_name="Sheet1", header=None) + +# ------------------------------------------------------------- +# READ STRUCTURE +# ------------------------------------------------------------- +scenario_row = df.iloc[0, 1:].tolist() +product_row = df.iloc[1, 1:].tolist() +components = df.iloc[2:, 0].astype(str).str.strip().tolist() +values = df.iloc[2:, 1:].astype(float) + +# ------------------------------------------------------------- +# CLEAN NANS +# ------------------------------------------------------------- +valid = [i for i, s in enumerate(scenario_row) if str(s) != "nan"] +scenario_row = [scenario_row[i] for i in valid] +product_row = [product_row[i] for i in valid] +values = values.iloc[:, valid] + +# ------------------------------------------------------------- +# MULTILINE SCENARIO LABELS (automatic wrapping) +# ------------------------------------------------------------- +scenario_row_wrapped = ["\n".join(textwrap.wrap(s, width=18)) for s in scenario_row] + +# ------------------------------------------------------------- +# BUILD MULTIINDEX +# ------------------------------------------------------------- +tuples = list(zip(scenario_row_wrapped, product_row)) +df_plot = pd.DataFrame(values.values, index=components, columns=pd.MultiIndex.from_tuples(tuples)) + +# ------------------------------------------------------------- +# FLATTENED PRODUCT LABELS +# ------------------------------------------------------------- +flat_products = product_row + +# ------------------------------------------------------------- +# GROUP POSITIONS FOR SCENARIO LABELS +# ------------------------------------------------------------- +scenario_groups = {} +for idx, scen in enumerate(scenario_row_wrapped): + scenario_groups.setdefault(scen, []).append(idx) + +x = np.arange(len(flat_products)) + +# ------------------------------------------------------------- +# COLOR MAP +# ------------------------------------------------------------- +color_map = { + "CapEx ($/kg)": "navy", + "OpEx ($/kg)": "orange", + "Feedstock ($/kg)": "deepskyblue", + "Taxes ($/kg)": "lightpink", + "Finances ($/kg)": "yellowgreen", +} + +# ------------------------------------------------------------- +# PLOT +# ------------------------------------------------------------- +plt.figure(figsize=(18, 7)) + +bottom = np.zeros(len(x)) + +for comp in components: + y = df_plot.loc[comp].values + plt.bar(x, y, bottom=bottom, color=color_map[comp], label=comp) + bottom += y + +# ------------------------------------------------------------- +# X-AXIS LABELS (PRODUCT LEVEL) +# ------------------------------------------------------------- +plt.xticks(x, flat_products, rotation=0, ha="center") + +# ------------------------------------------------------------- +# Y-AXIS LABEL +# ------------------------------------------------------------- +plt.ylabel("Levelized cost ($/kg)") + +plt.title("Cost Breakdown by Product and Scenario") + +# ------------------------------------------------------------- +# SCENARIO LABELS (CENTERED ABOVE GROUPS) +# ------------------------------------------------------------- +ymin, ymax = plt.ylim() +for scen, idxs in scenario_groups.items(): + center = np.mean(idxs) + plt.text( + center, ymax + ymax * 0.04, scen, ha="center", va="bottom", fontsize=11, fontweight="bold" + ) + +plt.ylim(ymin, ymax * 1.25) + +plt.legend(title="Cost Component", bbox_to_anchor=(1.02, 1), loc="upper left") +plt.tight_layout() +plt.show() diff --git a/examples/36_paper_mill/driver_config.yaml b/examples/36_paper_mill/driver_config.yaml new file mode 100644 index 000000000..bbba9ce27 --- /dev/null +++ b/examples/36_paper_mill/driver_config.yaml @@ -0,0 +1,10 @@ +name: driver_config +description: This analysis runs a paper mill plant and matches other examples in H2Integrate +general: + folder_output: outputs +recorder: + file: cases.sql + overwrite_recorder: true + flag: true + includes: ['*'] + excludes: ['*_resource*'] diff --git a/examples/36_paper_mill/paper_mill_profast_finance.py b/examples/36_paper_mill/paper_mill_profast_finance.py new file mode 100644 index 000000000..3d2a3b6de --- /dev/null +++ b/examples/36_paper_mill/paper_mill_profast_finance.py @@ -0,0 +1,198 @@ +import pytest +import openmdao.api as om +from pytest import fixture + +from h2integrate.finances.profast_lco import ProFastLCO + + +@fixture +def profast_inputs_no1(): + params = { + "analysis_start_year": 2030, # changed + "installation_time": 24, # 24 months? + "inflation_rate": 0.0, + "discount_rate": 0.0948, + "debt_equity_ratio": 1.72, + "property_tax_and_insurance": 0.015, # should we use this value? + "total_income_tax_rate": 0.2574, # should we use this value? + "capital_gains_tax_rate": 0.15, # should we use this value? + "sales_tax_rate": 0.00, + "debt_interest_rate": 0.046, # should we use this value? + "debt_type": "Revolving debt", + "loan_period_if_used": 0, + "cash_onhand_months": 1, + "admin_expense": 0.00, + } + cap_items = {"depr_type": "MACRS", "depr_period": 7, "refurb": [0.0]} + model_inputs = {"params": params, "capital_items": cap_items} + + return model_inputs + + +@fixture +def fake_filtered_tech_config(): + tech_config = { + "wind": {"model_inputs": {}}, + "solar": {"model_inputs": {}}, + "battery": {"model_inputs": {}}, + "natural_gas": {"model_inputs": {}}, + } + return tech_config + + +@fixture +def fake_cost_dict(): + fake_costs = { + "capex_adjusted_wind": 0, + "opex_adjusted_wind": 0, + "varopex_adjusted_wind": [0.0] * 0, + "capex_adjusted_solar": 0, + "opex_adjusted_solar": 0, + "varopex_adjusted_solar": [0.0] * 0, + "capex_adjusted_battery": 0, + "opex_adjusted_battery": 0, + "varopex_adjusted_battery": [0.0] * 00, + "capex_adjusted_natural_gas": 0, + "opex_adjusted_natural_gas": 0, + "varopex_adjusted_natural_gas": [0] * 00, + } + return fake_costs + + +@pytest.mark.regression +def test_profast_comp(profast_inputs_no1, fake_filtered_tech_config, fake_cost_dict, subtests): + mean_hourly_production = 34246.6 # ton/hr + prob = om.Problem() + plant_config = { + "plant": { + "plant_life": 40, + }, + "finance_parameters": {"model_inputs": profast_inputs_no1}, + } + pf = ProFastLCO( + driver_config={}, + plant_config=plant_config, + tech_config=fake_filtered_tech_config, + commodity_type="electricity", + description="no1", + ) + ivc = om.IndepVarComp() + + ivc.add_output("rated_electricity_production", mean_hourly_production, units="kW") + ivc.add_output("capacity_factor", [0.9] * plant_config["plant"]["plant_life"], units="unitless") + + prob.model.add_subsystem("ivc", ivc, promotes=["*"]) + prob.model.add_subsystem("pf", pf, promotes=["rated_electricity_production", "capacity_factor"]) + prob.setup() + for variable, cost in fake_cost_dict.items(): + units = "USD" if "capex" in variable else "USD/year" + prob.set_val(f"pf.{variable}", cost, units=units) + + prob.run_model() + + lcoe = prob.get_val("pf.LCOE_no1", units="USD/(MW*h)") + price = prob.get_val("pf.price_electricity_no1", units="USD/(MW*h)") + + wacc = prob.get_val("pf.wacc_electricity_no1", units="percent") + crf = prob.get_val("pf.crf_electricity_no1", units="percent") + profit_index = prob.get_val("pf.profit_index_electricity_no1", units="unitless") + irr = prob.get_val("pf.irr_electricity_no1", units="percent") + ipp = prob.get_val("pf.investor_payback_period_electricity_no1", units="yr") + + lcoe_breakdown = prob.get_val("pf.LCOE_no1_breakdown") + + with subtests.test("LCOE"): + assert pytest.approx(lcoe[0], rel=1e-6) == 63.8181779 + + with subtests.test("WACC"): + assert pytest.approx(wacc[0], rel=1e-6) == 0.056453864 + + with subtests.test("CRF"): + assert pytest.approx(crf[0], rel=1e-6) == 0.0674704169 + + with subtests.test("Profit Index"): + assert pytest.approx(profit_index[0], rel=1e-6) == 2.12026237778 + + with subtests.test("IRR"): + assert pytest.approx(irr[0], rel=1e-6) == 0.0948 + + with subtests.test("Investor payback period"): + assert pytest.approx(ipp[0], rel=1e-6) == 8 + + with subtests.test("LCOE == price"): + assert pytest.approx(lcoe, rel=1e-6) == price + + with subtests.test("LCOE breakdown total"): + assert pytest.approx(lcoe_breakdown["LCOE: Total ($/kWh)"] * 1e3, rel=1e-6) == lcoe + + +@pytest.mark.regression +def test_profast_comp_coproduct( + profast_inputs_no1, fake_filtered_tech_config, fake_cost_dict, subtests +): + mean_hourly_production = 500000.0 # kW*h + grid_sell_price = 63.8181779 / 1e3 # USD/(kW*h) + wind_sold_USD = [-1 * mean_hourly_production * 8760 * grid_sell_price] * 30 + fake_cost_dict.update({"varopex_adjusted_wind": wind_sold_USD}) + + prob = om.Problem() + plant_config = { + "plant": { + "plant_life": 30, + }, + "finance_parameters": {"model_inputs": profast_inputs_no1}, + } + pf = ProFastLCO( + driver_config={}, + plant_config=plant_config, + tech_config=fake_filtered_tech_config, + commodity_type="electricity", + description="no1", + ) + ivc = om.IndepVarComp() + ivc.add_output("rated_electricity_production", mean_hourly_production, units="kW") + ivc.add_output("capacity_factor", [1.0] * plant_config["plant"]["plant_life"], units="unitless") + + prob.model.add_subsystem("ivc", ivc, promotes=["*"]) + prob.model.add_subsystem("pf", pf, promotes=["rated_electricity_production", "capacity_factor"]) + prob.setup() + for variable, cost in fake_cost_dict.items(): + units = "USD" if "capex" in variable else "USD/year" + prob.set_val(f"pf.{variable}", cost, units=units) + + prob.run_model() + + lcoe = prob.get_val("pf.LCOE_no1", units="USD/(MW*h)") + price = prob.get_val("pf.price_electricity_no1", units="USD/(MW*h)") + + wacc = prob.get_val("pf.wacc_electricity_no1", units="percent") + crf = prob.get_val("pf.crf_electricity_no1", units="percent") + profit_index = prob.get_val("pf.profit_index_electricity_no1", units="unitless") + irr = prob.get_val("pf.irr_electricity_no1", units="percent") + ipp = prob.get_val("pf.investor_payback_period_electricity_no1", units="yr") + + lcoe_breakdown = prob.get_val("pf.LCOE_no1_breakdown") + + with subtests.test("LCOE"): + assert pytest.approx(lcoe[0], abs=1e-6) == 0 + + with subtests.test("WACC"): + assert pytest.approx(wacc[0], rel=1e-6) == 0.056453864 + + with subtests.test("CRF"): + assert pytest.approx(crf[0], rel=1e-6) == 0.0674704169 + + with subtests.test("Profit Index"): + assert pytest.approx(profit_index[0], rel=1e-6) == 2.12026237778 + + with subtests.test("IRR"): + assert pytest.approx(irr[0], rel=1e-6) == 0.0948 + + with subtests.test("Investor payback period"): + assert pytest.approx(ipp[0], rel=1e-6) == 8 + + with subtests.test("LCOE == price"): + assert pytest.approx(lcoe, rel=1e-6) == price + + with subtests.test("LCOE breakdown total"): + assert pytest.approx(lcoe_breakdown["LCOE: Total ($/kWh)"] * 1e3, rel=1e-6) == lcoe diff --git a/examples/36_paper_mill/plant_config.yaml b/examples/36_paper_mill/plant_config.yaml new file mode 100644 index 000000000..238882956 --- /dev/null +++ b/examples/36_paper_mill/plant_config.yaml @@ -0,0 +1,57 @@ +name: plant_config +description: This plant is located in MN, USA... +sites: + site: + latitude: 47.5233 + longitude: -92.5366 +# 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 +# add lignin as commodity +# paper, saf, lignin +technology_interconnections: + - [paper_mill, saf, lignin, pipe] +plant: + plant_life: 40 + simulation: + timezone: -6 +finance_parameters: + finance_groups: + finance_model: ProFastLCO + model_inputs: + params: + analysis_start_year: 2030 + installation_time: 24 # months + inflation_rate: 0.0 # 0 for nominal analysis + discount_rate: 0.09 # nominal return based on 2024 ATB baseline workbook for land-based wind + debt_equity_ratio: 1.86 # 2024 ATB uses 72.4% debt for land-based wind + property_tax_and_insurance: 0.03 # percent of CAPEX estimated based on https://www.nrel.gov/docs/fy25osti/91775.pdf https://www.house.mn.gov/hrd/issinfo/clsrates.aspx + total_income_tax_rate: 0.308 # 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.07375 # total state and local sales tax in St. Louis County https://taxmaps.state.mn.us/salestax/ + debt_interest_rate: 0.06 # 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.] + finance_subgroups: + paper: + commodity: paper + commodity_stream: paper_mill + technologies: [paper_mill] + pulp: + commodity: pulp_out + commodity_stream: paper_mill + technologies: [paper_mill] + saf: + commodity: saf + commodity_stream: saf + technologies: [saf] + cost_adjustment_parameters: + cost_year_adjustment_inflation: 0.0 # used to adjust modeled costs to target_dollar_year + target_dollar_year: 2022 diff --git a/examples/36_paper_mill/run_paper_mill_mn.py b/examples/36_paper_mill/run_paper_mill_mn.py new file mode 100644 index 000000000..a76122491 --- /dev/null +++ b/examples/36_paper_mill/run_paper_mill_mn.py @@ -0,0 +1,31 @@ +"""Run script for the paper mill example (Minnesota site). + +Sets up, runs, and post-processes the H2Integrate paper mill model, then exports +results to CSV using the built-in SQL postprocessing utilities: + +- Scalar outputs are saved via ``convert_sql_to_csv_summary``. +- Timeseries profiles are saved via ``save_case_timeseries_as_csv``. + +Both CSV files are written to the ``outputs/`` directory alongside the SQL recorder file. +""" + +import os + +from h2integrate import EXAMPLE_DIR +from h2integrate.core.h2integrate_model import H2IntegrateModel +from h2integrate.postprocess.sql_to_csv import convert_sql_to_csv_summary +from h2integrate.postprocess.sql_timeseries_to_csv import save_case_timeseries_as_csv + + +os.chdir(EXAMPLE_DIR / "36_paper_mill") + +model = H2IntegrateModel("36_paper_mill_mn.yaml") + +model.setup() +model.run() +model.post_process() + +sql_fpath = EXAMPLE_DIR / "36_paper_mill" / "outputs" / "cases.sql" + +convert_sql_to_csv_summary(sql_fpath) +save_case_timeseries_as_csv(sql_fpath) diff --git a/examples/36_paper_mill/tech_config.yaml b/examples/36_paper_mill/tech_config.yaml new file mode 100644 index 000000000..87ed3acd7 --- /dev/null +++ b/examples/36_paper_mill/tech_config.yaml @@ -0,0 +1,75 @@ +# TOASK: (1) Easy way to switch on and off SAF portion? Create another example?; (2) paper vs pulp +name: technology_config +description: This plant produces pulp and paper +technologies: + paper_mill: + performance_model: + model: PaperMillPerformanceModel + cost_model: + model: PaperMillCostModel + model_inputs: + shared_parameters: + plant_capacity_mtpy: 3300000 + performance_parameters: + capacity_factor: 0.9 + cost_parameters: + operational_year: 2030 + installation_time: 24 # months + inflation_rate: 0.0 # 0 for nominal analysis + # Feedstock parameters (flattened) + wood_unitcost: 99.3 + wood_transport_cost: 0.0 + electricity_cost: 0.054 + calcium_carbonate_unitcost: 0.33 #$/kg consumable + sodium_sulfide_unitcost: 0.22 #$/kg consumable + sodium_hydroxide_unitcost: 0.33 #$/kg consumable + chlorine_dioxide_unitcost: 1.79 #$/kg consumable + hydrogen_peroxide_unitcost: 0.33 #$/kg consumable + magnesium_sulfate_unitcost: 0.36 #$/kg consumable + oxygen_unitcost: 0 #$/ton consumable + raw_water_unitcost: 0.001519 + raw_water_consumption: 40693 + wood_consumption: 0.225 + calcium_carbonate_consumption: 240 # kg/MT product + sodium_sulfide_consumption: 16.5 # kg/MT product + sodium_hydroxide_consumption: 3.5 # kg/MT product + chlorine_dioxide_consumption: 2.0 # kg/MT product + hydrogen_peroxide_consumption: .5 # kg/MT product + magnesium_sulfate_consumption: .2 # kg/MT product + oxygen_consumption: 2.25 # kg/MT product + electricity_consumption: 68.7 + water_disposal_unitcost: 0.002013 + water_disposal_rate: 18927 + saf: + performance_model: + model: SAFPerformanceModel + cost_model: + model: SAFCostModel + model_inputs: + shared_parameters: + plant_capacity_mtpy: 28000 + performance_parameters: + capacity_factor: 0.9 + cost_parameters: + operational_year: 2030 +# TOASK: LCOH vs hydrogen unit cost. Should we use the hydrogen module for the August results +# lcoh: 7.37 + installation_time: 24 # months + inflation_rate: 0.0 # 0 for nominal analysis + # Feedstock parameters (flattened) + lignin_unitcost: 0.78 + lignin_transport_cost: 0.0 + lignin_consumption: 1650 + hydrogen_unitcost: 7.37 + hydrogen_transport_cost: 0.0 + hydrogen_consumption: 580 + electricity_cost: 0.054 + electricity_consumption: 19750 + raw_water_unitcost: 0.001519 + raw_water_consumption: 2839 + water_disposal_unitcost: 0.002013 + water_disposal_rate: 0 + salt_mix_unitcost: 0.86 # $/kg + salt_mix_consumption: 41.3 # kg/t SAF + hydrogen_chloride_unitcost: 0.26 # $/kg + hydrogen_chloride_consumption: 1.5 # kg/t SAF diff --git a/h2integrate/converters/paper_mill/__init__.py b/h2integrate/converters/paper_mill/__init__.py new file mode 100644 index 000000000..0ff4f9f6b --- /dev/null +++ b/h2integrate/converters/paper_mill/__init__.py @@ -0,0 +1,4 @@ +from h2integrate.converters.paper_mill.paper_mill import ( + PaperMillPerformanceModel, + PaperMillCostModel, +) diff --git a/h2integrate/converters/paper_mill/paper_mill.py b/h2integrate/converters/paper_mill/paper_mill.py new file mode 100644 index 000000000..64f946fa1 --- /dev/null +++ b/h2integrate/converters/paper_mill/paper_mill.py @@ -0,0 +1,201 @@ +from attrs import field, define + +from h2integrate.core.utilities import BaseConfig, merge_shared_inputs +from h2integrate.core.validators import must_equal +from h2integrate.core.model_baseclasses import CostModelBaseClass, PerformanceModelBaseClass + + +@define(kw_only=True) +class PaperMillPerformanceModelConfig(BaseConfig): + plant_capacity_mtpy: float = field() + capacity_factor: float = field() + + +class PaperMillPerformanceModel(PerformanceModelBaseClass): + """ + An OpenMDAO component for modeling the performance of an paper mill plant. + Computes annual paper production based on plant capacity and capacity factor. + """ + + _time_step_bounds = ( + 3600, + 3600, + ) # (min, max) time step lengths (in seconds) compatible with this model + + def initialize(self): + super().initialize() + self.commodity = "paper" + self.commodity_amount_units = "t" + self.commodity_rate_units = "t/h" + + def setup(self): + super().setup() + + self.config = PaperMillPerformanceModelConfig.from_dict( + merge_shared_inputs(self.options["tech_config"]["model_inputs"], "performance"), + additional_cls_name=self.__class__.__name__, + ) + + self.add_input("plant_capacity_mtpy", val=self.config.plant_capacity_mtpy, units="t/year") + + n_timesteps = self.options["plant_config"]["plant"]["simulation"]["n_timesteps"] + self.add_output("lignin_out", val=0.0, shape=n_timesteps, units="kg/h") + self.add_output("rated_lignin_production", shape=n_timesteps, val=0.0, units="kg/h") + self.add_output("total_lignin_produced", val=0.0, units="kg") + self.add_output("annual_lignin_produced", val=0.0, units="kg/year") + + self.add_output("pulp_out", val=0.0, shape=n_timesteps, units="t/h") + self.add_output("rated_pulp_out_production", val=0.0, units="t/h") + self.add_output("total_pulp_out_produced", val=0.0, units="t") + self.add_output("annual_pulp_out_produced", val=0.0, units="t/year") + + def compute(self, inputs, outputs): + plant_capacity_mtpy = inputs["plant_capacity_mtpy"] + capacity_factor = self.config.capacity_factor + paper_mill_production_mtpy = plant_capacity_mtpy * capacity_factor + outputs["paper_out"] = paper_mill_production_mtpy / 8760 # tons per hour + outputs["rated_paper_production"] = plant_capacity_mtpy / 8760 # tons per hour + outputs["capacity_factor"] = capacity_factor + outputs["total_paper_produced"] = outputs["paper_out"].sum() + outputs["annual_paper_produced"] = outputs["total_paper_produced"] * ( + 1 / self.fraction_of_year_simulated + ) + + outputs["lignin_out"] = 0.06 * 1000 * paper_mill_production_mtpy / 8760 # tons per hour + outputs["rated_lignin_production"] = 0.06 * 1000 * plant_capacity_mtpy / 8760 + outputs["total_lignin_produced"] = outputs["lignin_out"].sum() + outputs["annual_lignin_produced"] = outputs["total_lignin_produced"] * ( + 1 / self.fraction_of_year_simulated + ) + + outputs["pulp_out"] = 1.1 * paper_mill_production_mtpy / 8760 + outputs["rated_pulp_out_production"] = 1.1 * plant_capacity_mtpy / 8760 + outputs["total_pulp_out_produced"] = outputs["pulp_out"].sum() + outputs["annual_pulp_out_produced"] = outputs["total_pulp_out_produced"] * ( + 1 / self.fraction_of_year_simulated + ) + + +@define(kw_only=True) +class PaperMillCostModelConfig(BaseConfig): + installation_time: int = field() + inflation_rate: float = field() + operational_year: int = field() + plant_capacity_mtpy: float = field() + cost_year: int = field(default=2023, converter=int, validator=must_equal(2023)) + wood_unitcost: float = field(default=99.3) # $/MT of wood + wood_transport_cost: float = field(default=0.0) + + calcium_carbonate_unitcost: float = field(default=0.33) # $/kg consumable + calcium_carbonate_transport_cost: float = field(default=0.0) + sodium_sulfide_unitcost: float = field(default=0.22) # $/kg consumable + sodium_sulfide_transport_cost: float = field(default=0.0) + sodium_hydroxide_unitcost: float = field(default=0.33) # $/kg consumable + sodium_hydroxide_transport_cost: float = field(default=0.0) + chlorine_dioxide_unitcost: float = field(default=1.79) # $/kg consumable + chlorine_dioxide_transport_cost: float = field(default=0.0) + hydrogen_peroxide_unitcost: float = field(default=0.33) # $/kg consumable + hydrogen_peroxide_transport_cost: float = field(default=0.0) + magnesium_sulfate_unitcost: float = field(default=0.36) # $/kg consumable + magnesium_sulfate_transport_cost: float = field(default=0.0) + oxygen_unitcost: float = field(default=0) # $/ton consumable + oxygen_transport_cost: float = field(default=0.0) + + electricity_cost: float = field(default=0.054) # $/kWh + raw_water_unitcost: float = field(default=0.001519) # $/kg water + wood_consumption: float = field(default=0.225) # MT/MT product + raw_water_consumption: float = field(default=40693) # kg/tonne product + + calcium_carbonate_consumption: float = field(default=240) # kg/MT product + sodium_sulfide_consumption: float = field(default=16.5) # kg/MT product + sodium_hydroxide_consumption: float = field(default=3.5) # kg/MT product + chlorine_dioxide_consumption: float = field(default=2) # kg/MT product + hydrogen_peroxide_consumption: float = field(default=0.5) # kg/MT product + magnesium_sulfate_consumption: float = field(default=0.2) # kg/MT product + oxygen_consumption: float = field(default=2.25) # kg/MT product + + electricity_consumption: float = field(default=68.7) # kWh/tonne product + water_disposal_unitcost: float = field(default=0.002013) # $/kg + water_disposal_rate: float = field(default=18927) # kg/MT product + + +class PaperMillCostModel(CostModelBaseClass): + """ + An OpenMDAO component for calculating the costs associated with paper mill production. + Includes CapEx, OpEx, and byproduct credits. + """ + + _time_step_bounds = ( + 3600, + 3600, + ) # (min, max) time step lengths (in seconds) compatible with this model + + def setup(self): + self.config = PaperMillCostModelConfig.from_dict( + merge_shared_inputs(self.options["tech_config"]["model_inputs"], "cost"), + additional_cls_name=self.__class__.__name__, + ) + super().setup() + + self.add_input("plant_capacity_mtpy", val=self.config.plant_capacity_mtpy, units="t/year") + + def compute(self, inputs, outputs, discrete_inputs, discrete_outputs): + plant_capacity_mtpy = inputs["plant_capacity_mtpy"][0] + + # Calculate plant CapEx for Kraft process + total_plant_capex = 2500 * plant_capacity_mtpy + + # Fixed O&M Costs + # TODO: Need to update labor cost + labor_cost_annual_operation = ( + 69375996.9 + * ((plant_capacity_mtpy / 365 * 1000) ** 0.25242) + / ((1162077 / 365 * 1000) ** 0.25242) + ) + labor_cost_maintenance = 0.00863 * total_plant_capex + 0.25 * (labor_cost_annual_operation + labor_cost_maintenance) + + fixed_operating_cost = 370 * plant_capacity_mtpy + + property_tax_insurance = 0.02 * total_plant_capex + + total_fixed_operating_cost = fixed_operating_cost + property_tax_insurance + + # Sum total consumables costs and consumption + c = self.config + consumable_costs_per_mt = { + "raw_water": c.raw_water_consumption * c.raw_water_unitcost, + "wood": c.wood_consumption * (c.wood_unitcost + c.wood_transport_cost), + "calcium_carbonate": c.calcium_carbonate_consumption + * (c.calcium_carbonate_unitcost + c.calcium_carbonate_transport_cost), + "sodium_sulfide": c.sodium_sulfide_consumption + * (c.sodium_sulfide_unitcost + c.sodium_sulfide_transport_cost), + "sodium_hydroxide": c.sodium_hydroxide_consumption + * (c.sodium_hydroxide_unitcost + c.sodium_hydroxide_transport_cost), + "chlorine_dioxide": c.chlorine_dioxide_consumption + * (c.chlorine_dioxide_unitcost + c.chlorine_dioxide_transport_cost), + "hydrogen_peroxide": c.hydrogen_peroxide_consumption + * (c.hydrogen_peroxide_unitcost + c.hydrogen_peroxide_transport_cost), + "magnesium_sulfate": c.magnesium_sulfate_consumption + * (c.magnesium_sulfate_unitcost + c.magnesium_sulfate_transport_cost), + "oxygen": c.oxygen_consumption * (c.oxygen_unitcost + c.oxygen_transport_cost), + } + variable_consumables_cost = plant_capacity_mtpy * sum(consumable_costs_per_mt.values()) + + water_disposal_cost = ( + plant_capacity_mtpy + * self.config.water_disposal_unitcost + * self.config.water_disposal_rate + ) + + electricity_cost = plant_capacity_mtpy * ( + self.config.electricity_consumption * self.config.electricity_cost + ) + + total_variable_operating_cost = ( + variable_consumables_cost + water_disposal_cost + electricity_cost + ) + + outputs["CapEx"] = total_plant_capex + outputs["OpEx"] = total_fixed_operating_cost + outputs["VarOpEx"] = total_variable_operating_cost diff --git a/h2integrate/converters/saf/__init__.py b/h2integrate/converters/saf/__init__.py new file mode 100644 index 000000000..09202ea7f --- /dev/null +++ b/h2integrate/converters/saf/__init__.py @@ -0,0 +1,4 @@ +from h2integrate.converters.saf.saf import ( + SAFPerformanceModel, + SAFCostModel, +) diff --git a/h2integrate/converters/saf/saf.py b/h2integrate/converters/saf/saf.py new file mode 100644 index 000000000..bbd87f022 --- /dev/null +++ b/h2integrate/converters/saf/saf.py @@ -0,0 +1,150 @@ +from attrs import field, define + +from h2integrate.core.utilities import BaseConfig, merge_shared_inputs +from h2integrate.core.validators import must_equal +from h2integrate.core.model_baseclasses import CostModelBaseClass, PerformanceModelBaseClass + + +@define(kw_only=True) +class SAFPerformanceModelConfig(BaseConfig): + plant_capacity_mtpy: float = field() + capacity_factor: float = field() + + +class SAFPerformanceModel(PerformanceModelBaseClass): + """ + An OpenMDAO component for modeling the performance of a saf plant. + Computes annual saf production based on plant capacity and capacity factor. + """ + + _time_step_bounds = ( + 3600, + 3600, + ) # (min, max) time step lengths (in seconds) compatible with this model + + def initialize(self): + super().initialize() + self.commodity = "saf" + self.commodity_amount_units = "t" + self.commodity_rate_units = "t/h" + + def setup(self): + super().setup() + self.config = SAFPerformanceModelConfig.from_dict( + merge_shared_inputs(self.options["tech_config"]["model_inputs"], "performance"), + additional_cls_name=self.__class__.__name__, + ) + + self.add_input("plant_capacity_mtpy", val=self.config.plant_capacity_mtpy, units="t/year") + n_timesteps = self.options["plant_config"]["plant"]["simulation"]["n_timesteps"] + self.add_input("lignin_in", val=0.0, shape=n_timesteps, units="kg/h") + + def compute(self, inputs, outputs): + plant_capacity_mtpy = inputs["plant_capacity_mtpy"] + capacity_factor = self.config.capacity_factor + saf_production_mtpy = plant_capacity_mtpy * capacity_factor + outputs["saf_out"] = saf_production_mtpy / 8760 + outputs["rated_saf_production"] = plant_capacity_mtpy / 8760 + outputs["capacity_factor"] = capacity_factor + outputs["total_saf_produced"] = outputs["saf_out"].sum() + outputs["annual_saf_produced"] = outputs["total_saf_produced"] * ( + 1 / self.fraction_of_year_simulated + ) + + +@define(kw_only=True) +class SAFCostModelConfig(BaseConfig): + installation_time: int = field() + inflation_rate: float = field() + operational_year: int = field() + plant_capacity_mtpy: float = field() + cost_year: int = field(default=2023, converter=int, validator=must_equal(2023)) + + # Feedstock parameters - flattened from the nested structure + lignin_unitcost: float = field(default=0.78) # $/kg of final product + lignin_transport_cost: float = field(default=0.0) + salt_mix_unitcost: float = field(default=0.86) # $/kg consumable + salt_mix_transport_cost: float = field(default=0.0) + hydrogen_chloride_unitcost: float = field(default=0.26) # $/kg consumable + hydrogen_chloride_transport_cost: float = field(default=0.0) + hydrogen_unitcost: float = field(default=7.37) # $/kg consumable + hydrogen_transport_cost: float = field(default=0.0) + electricity_cost: float = field(default=0.054) # $/kWh + raw_water_unitcost: float = field(default=0.001519) # $/kg water + lignin_consumption: float = field(default=1650) # kg/MT product + raw_water_consumption: float = field(default=2839) # kg/tonne product + hydrogen_consumption: float = field(default=580) # kg/tonne product + salt_mix_consumption: float = field(default=41.3) # kg/MT product + hydrogen_chloride_consumption: float = field(default=1.5) # kg/MT product + electricity_consumption: float = field(default=19750) # kWh/tonne product + water_disposal_unitcost: float = field(default=0.002013) # $/kg + water_disposal_rate: float = field(default=0) # TODO: Change assumption + + +class SAFCostModel(CostModelBaseClass): + """ + An OpenMDAO component for calculating the costs associated with saf production. + Includes CapEx, OpEx, and byproduct credits. + """ + + _time_step_bounds = ( + 3600, + 3600, + ) # (min, max) time step lengths (in seconds) compatible with this model + + # TOASK: In that case, do we need this function? + def setup(self): + self.config = SAFCostModelConfig.from_dict( + merge_shared_inputs(self.options["tech_config"]["model_inputs"], "cost"), + additional_cls_name=self.__class__.__name__, + ) + super().setup() + + self.add_input("plant_capacity_mtpy", val=self.config.plant_capacity_mtpy, units="t/year") + + def compute(self, inputs, outputs, discrete_inputs, discrete_outputs): + plant_capacity_mtpy = inputs["plant_capacity_mtpy"][0] + + # Calculate saf production costs directly + total_plant_capex = 5570 * plant_capacity_mtpy + + # Fixed O&M Costs + # TODO: Need to update labor cost + labor_cost_annual_operation = ( + 69375996.9 + * ((plant_capacity_mtpy / 365 * 1000) ** 0.25242) + / ((1162077 / 365 * 1000) ** 0.25242) + ) + labor_cost_maintenance = 0.00863 * total_plant_capex + 0.25 * (labor_cost_annual_operation + labor_cost_maintenance) + + fixed_operating_cost = 390 * plant_capacity_mtpy + + property_tax_insurance = 0.02 * total_plant_capex + + total_fixed_operating_cost = fixed_operating_cost + property_tax_insurance + + c = self.config + consumable_costs_per_mt = { + "raw_water": c.raw_water_consumption * c.raw_water_unitcost, + "lignin": c.lignin_consumption * (c.lignin_unitcost + c.lignin_transport_cost), + "salt_mix": c.salt_mix_consumption * (c.salt_mix_unitcost + c.salt_mix_transport_cost), + "hydrogen_chloride": c.hydrogen_chloride_consumption + * (c.hydrogen_chloride_unitcost + c.hydrogen_chloride_transport_cost), + "hydrogen": c.hydrogen_consumption * (c.hydrogen_unitcost + c.hydrogen_transport_cost), + } + variable_consumables_cost = plant_capacity_mtpy * sum(consumable_costs_per_mt.values()) + + water_disposal_cost = ( + plant_capacity_mtpy * c.water_disposal_unitcost * c.water_disposal_rate + ) + + electricity_cost = plant_capacity_mtpy * (c.electricity_consumption * c.electricity_cost) + + total_variable_operating_cost = ( + variable_consumables_cost + water_disposal_cost + electricity_cost + ) + + outputs["CapEx"] = total_plant_capex + outputs["OpEx"] = total_fixed_operating_cost + outputs["VarOpEx"] = total_variable_operating_cost diff --git a/h2integrate/core/supported_models.py b/h2integrate/core/supported_models.py index 51fa624bd..f98015a29 100644 --- a/h2integrate/core/supported_models.py +++ b/h2integrate/core/supported_models.py @@ -133,6 +133,10 @@ def copy(self): "SimpleThermalNuclearReactorCostModel": "converters.nuclear:SimpleThermalNuclearReactorCostModel", "SimpleThermalNuclearReactorPerformanceModel": "converters.nuclear:SimpleThermalNuclearReactorPerformanceModel", "NaturalGasCostModel": "converters.natural_gas:NaturalGasCostModel", + "PaperMillPerformanceModel": "converters.paper_mill:PaperMillPerformanceModel", + "PaperMillCostModel": "converters.paper_mill:PaperMillCostModel", + "SAFPerformanceModel": "converters.saf:SAFPerformanceModel", + "SAFCostModel": "converters.saf:SAFCostModel", # Transport "cable": "transporters:CablePerformanceModel", "pipe": "transporters:PipePerformanceModel", diff --git a/h2integrate/transporters/pipe.py b/h2integrate/transporters/pipe.py index 5ad3a92c1..60f2e22ff 100644 --- a/h2integrate/transporters/pipe.py +++ b/h2integrate/transporters/pipe.py @@ -24,6 +24,7 @@ def initialize(self): "wellhead_gas", "water", "oxygen", + "lignin", ], ) self.options.declare("plant_config", types=dict)