Merge pull request #312 from jaredthomas68/uniform-output

Uniform output

greenheart/simulation/greenheart_simulation.py

@@ -8,11 +8,20 @@
 
 pd.options.mode.chained_assignment = None  # default='warn'
 
+from hopp.simulation import HoppInterface
+from ProFAST import ProFAST
+
 from greenheart.simulation.technologies.ammonia.ammonia import (
-    run_ammonia_full_model,
+    run_ammonia_full_model, 
+    AmmoniaCostModelOutputs,
+    AmmoniaFinanceModelOutputs, 
+    AmmoniaCapacityModelOutputs
 )
 from greenheart.simulation.technologies.steel.steel import (
     run_steel_full_model,
+    SteelCostModelOutputs,
+    SteelFinanceModelOutputs,
+    SteelCapacityModelOutputs
 )
 
 # visualization imports
@@ -161,6 +170,66 @@ def __attrs_post_init__(self):
                     self.orbit_config["plant"]["capacity"] * 1e6
                 )
 
+@define
+class GreenHeartSimulationOutput:
+    """This is a dataclass to contain the outputs from GreenHEART
+
+        Args:
+            greenheart_config (GreenHeartSimulationConfig): all inputs to the greenheart simulation
+            hopp_interface (HoppInterface): the hopp interface created and used by GreenHEART in the simulation
+            profast_lcoe (ProFAST): the profast instance used for the lcoe calculations
+            profast_lcoh (ProFAST): the profast instance used for the lcoh calculations
+            profast_lcoh (ProFAST): the profast instance used for the lcoh calculations if  hydrogen were produced only from the grid
+            lcoe (float): levelized cost of energy (electricity)
+            lcoh (float): levelized cost of hydrogen
+            lcoh_grid_only (float): levelized cost of hydrogen if produced only from the grid
+            hopp_results (dict): results from the hopp simulation
+            electrolyzer_physics_results (dict): results of the electrolysis simulation
+            capex_breakdown (dict): overnight capex broken down by technology
+            opex_breakdown_annual (dict): annual operational expenditures broken down by technology
+            annual_energy_breakdown (dict): annual energy generation and usage broken down by technology
+            hourly_energy_breakdown (dict): hourly energy generation and usage broken down by technology
+            remaining_power_profile (np.ndarray): unused power (hourly)
+            steel_capacity (Optional[SteelCapacityModelOutputs]): steel capacity information
+            steel_costs (Optional[SteelCostModelOutputs]): steel cost information
+            steel_finance (Optional[SteelFinanceModelOutputs]): steel financial information
+            ammonia_capacity (Optional[AmmoniaCapacityModelOutputs]): ammonia capacity information
+            ammonia_costs (Optional[AmmoniaCostModelOutputs]): ammonia cost information
+            ammonia_finance (Optional[AmmoniaFinanceModelOutputs]): ammonia finance information
+    """
+
+    # detailed simulation information
+    greenheart_config: GreenHeartSimulationConfig
+    hopp_interface: HoppInterface
+
+    # detailed financial outputs
+    profast_lcoe: ProFAST
+    profast_lcoh: ProFAST
+    profast_lcoh_grid_only: ProFAST
+
+    # high-level results
+    lcoe: float
+    lcoh: float
+    lcoh_grid_only: float
+
+    # detailed output information
+    hopp_results: dict
+    electrolyzer_physics_results: dict
+    capex_breakdown: dict
+    opex_breakdown_annual: dict
+    annual_energy_breakdown: dict
+    hourly_energy_breakdown: dict
+    remaining_power_profile: np.ndarray
+
+    # optional outputs
+    steel_capacity: Optional[SteelCapacityModelOutputs] = field(default=None)
+    steel_costs: Optional[SteelCostModelOutputs] = field(default=None)
+    steel_finance: Optional[SteelFinanceModelOutputs] = field(default=None)
+
+    ammonia_capacity: Optional[AmmoniaCapacityModelOutputs] = field(default=None)
+    ammonia_costs: Optional[AmmoniaCostModelOutputs] = field(default=None)
+    ammonia_finance: Optional[AmmoniaFinanceModelOutputs] = field(default=None)
+
 def setup_greenheart_simulation(config: GreenHeartSimulationConfig):
 
     # run orbit for wind plant construction and other costs
@@ -520,7 +589,7 @@ def energy_internals(
                 h2_storage_results,
                 total_accessory_power_renewable_kw,
                 total_accessory_power_grid_kw,
-                remaining_power_profile,
+                remaining_power_profile
             )
 
     # define function to provide to the brent solver
@@ -733,7 +802,7 @@ def simple_solver(initial_guess=0.0):
                     "hydrogen_amount_kgpy"
                 ] = hydrogen_amount_kgpy
 
-            _, _, steel_finance = run_steel_full_model(config.greenheart_config, save_plots=config.save_plots, show_plots=config.show_plots, output_dir=config.output_dir, design_scenario_id=config.design_scenario["id"])
+            steel_capacity, steel_costs, steel_finance = run_steel_full_model(config.greenheart_config, save_plots=config.save_plots, show_plots=config.show_plots, output_dir=config.output_dir, design_scenario_id=config.design_scenario["id"])
 
         else:
             steel_finance = {}
@@ -751,14 +820,14 @@ def simple_solver(initial_guess=0.0):
                     "hydrogen_amount_kgpy"
                 ] = hydrogen_amount_kgpy
 
-            _, _, ammonia_finance = run_ammonia_full_model(config.greenheart_config, save_plots=config.save_plots, show_plots=config.show_plots, output_dir=config.output_dir, design_scenario_id=config.design_scenario["id"])
+            ammonia_capacity, ammonia_costs, ammonia_finance = run_ammonia_full_model(config.greenheart_config, save_plots=config.save_plots, show_plots=config.show_plots, output_dir=config.output_dir, design_scenario_id=config.design_scenario["id"])
 
         else:
             ammonia_finance = {}
 
     ################# end OSW intermediate calculations
     if config.post_processing:
-        power_breakdown = he_util.post_process_simulation(
+        annual_energy_breakdown, hourly_energy_breakdown = he_util.post_process_simulation(
             lcoe,
             lcoh,
             pf_lcoh,
@@ -808,7 +877,7 @@ def simple_solver(initial_guess=0.0):
             pf_lcoh,
             electrolyzer_physics_results,
             pf_lcoe,
-            power_breakdown,
+            annual_energy_breakdown,
         )
     elif config.output_level == 4:
         return lcoe, lcoh, lcoh_grid_only
@@ -818,7 +887,30 @@ def simple_solver(initial_guess=0.0):
         return hopp_results, electrolyzer_physics_results, remaining_power_profile
     elif config.output_level == 7:
         return lcoe, lcoh, steel_finance, ammonia_finance
-
+    elif config.output_level == 8:
+        return GreenHeartSimulationOutput(
+            config,
+            hi,
+            pf_lcoe,
+            pf_lcoh,
+            pf_grid_only,
+            lcoe,
+            lcoh,
+            lcoh_grid_only,
+            hopp_results,
+            electrolyzer_physics_results,
+            capex_breakdown,
+            opex_breakdown_annual,
+            annual_energy_breakdown,
+            hourly_energy_breakdown,
+            remaining_power_profile,
+            steel_capacity = None if "steel" not in config.greenheart_config else steel_capacity, 
+            steel_costs = None if "steel" not in config.greenheart_config else steel_costs, 
+            steel_finance = None if "steel" not in config.greenheart_config else steel_finance,
+            ammonia_capacity = None if "ammonia" not in config.greenheart_config else ammonia_capacity, 
+            ammonia_costs = None if "ammonia" not in config.greenheart_config else ammonia_costs, 
+            ammonia_finance = None if "ammonia" not in config.greenheart_config else ammonia_finance
+        )
 
 def run_sweeps(simulate=False, verbose=True, show_plots=True, use_profast=True, output_dir="output/"):
 

greenheart/tools/eco/electrolysis.py

@@ -190,7 +190,7 @@ def run_electrolyzer_physics(
         ax[0, 0].plot(wind_speed)
         convolved_wind_speed = np.convolve(wind_speed, np.ones(N) / (N), mode="valid")
         ave_x = range(N, len(convolved_wind_speed) + N)
-
+        
         ax[0, 1].plot(ave_x, convolved_wind_speed)
         ax[0, 0].set(ylabel="Wind\n(m/s)", ylim=[0, 30], xlim=[0, len(wind_speed)])
         tick_spacing = 10
@@ -199,11 +199,14 @@ def run_electrolyzer_physics(
         y = greenheart_config["electrolyzer"]["rating"]
         ax[1, 0].plot(energy_to_electrolyzer_kw * 1e-3)
         ax[1, 0].axhline(y=y, color="r", linestyle="--", label="Nameplate Capacity")
-        ax[1, 1].plot(
-            ave_x[:-1],
-            np.convolve(
+
+        convolved_energy_to_electrolyzer = np.convolve(
                 energy_to_electrolyzer_kw * 1e-3, np.ones(N) / (N), mode="valid"
-            ),
+            )
+
+        ax[1, 1].plot(
+            ave_x,
+            convolved_energy_to_electrolyzer,
         )
         ax[1, 1].axhline(y=y, color="r", linestyle="--", label="Nameplate Capacity")
         ax[1, 0].set(
@@ -220,16 +223,10 @@ def run_electrolyzer_physics(
             ]
             * 1e-3
         )
+        convolved_hydrogen_production = np.convolve(electrolyzer_physics_results["H2_Results"]["Hydrogen Hourly Production [kg/hr]"]*1e-3,np.ones(N) / (N), mode="valid")
         ax[2, 1].plot(
-            ave_x[:-1],
-            np.convolve(
-                electrolyzer_physics_results["H2_Results"][
-                    "Hydrogen Hourly Production [kg/hr]"
-                ]
-                * 1e-3,
-                np.ones(N) / (N),
-                mode="valid",
-            ),
+            ave_x,
+            convolved_hydrogen_production,
         )
         tick_spacing = 2
         ax[2, 0].set(

greenheart/tools/eco/hopp_mgmt.py

@@ -199,7 +199,7 @@ def run_hopp(hi, project_lifetime, verbose=False):
         "hopp_interface": hi,
         "hybrid_plant": hi.system,
         "combined_hybrid_power_production_hopp": \
-            hi.system.grid._system_model.Outputs.system_pre_interconnect_kwac[0:8759],
+            hi.system.grid._system_model.Outputs.system_pre_interconnect_kwac[0:8760],
         "combined_hybrid_curtailment_hopp": hi.system.grid.generation_curtailed,
         "energy_shortfall_hopp": hi.system.grid.missed_load,
         "annual_energies": hi.system.annual_energies,

greenheart/tools/eco/utilities.py

@@ -794,11 +794,11 @@ def add_turbines(ax, turbine_x, turbine_y, radius, color):
         e_side_x = electrolyzer_area / e_side_y
         d_side_y = equipment_platform_side_length
         d_side_x = desal_equipment_area / d_side_y
-        ex = dx + d_side_x
-        ey = dy
+        electrolyzer_x = dx + d_side_x
+        electrolyzer_y = dy
 
         electrolyzer_patch = patches.Rectangle(
-            (ex, ey),
+            (electrolyzer_x, electrolyzer_y),
             e_side_x,
             e_side_y,
             color=electrolyzer_color,
@@ -1293,10 +1293,12 @@ def add_turbines(ax, turbine_x, turbine_y, radius, color):
     return 0
 
 
-def save_power_series(
-    hybrid_plant: HoppInterface.system, ax=None, simulation_length=8760, output_dir="./output/",
+def save_energy_flows(
+    hybrid_plant: HoppInterface.system, electrolyzer_physics_results, solver_results, hours, ax=None, simulation_length=8760, output_dir="./output/",
 ):
 
+
+
     if ax == None:
         fig, ax = plt.subplots(1)
 
@@ -1305,30 +1307,41 @@ def save_power_series(
         solar_plant_power = np.array(
             hybrid_plant.pv.generation_profile[0:simulation_length]
         )
-        output.update({"pv": solar_plant_power})
+        output.update({"pv generation [kW]": solar_plant_power})
     if hybrid_plant.wind:
         wind_plant_power = np.array(
             hybrid_plant.wind.generation_profile[0:simulation_length]
         )
-        output.update({"wind": wind_plant_power})
+        output.update({"wind generation [kW]": wind_plant_power})
     if hybrid_plant.wave:
         wave_plant_power = np.array(
             hybrid_plant.wave.generation_profile[0:simulation_length]
         )
-        output.update({"wave": wave_plant_power})
+        output.update({"wave generation [kW]": wave_plant_power})
     if hybrid_plant.battery:
-        battery_power_out = hybrid_plant.battery.outputs.dispatch_P
-        output.update({"battery": battery_power_out})
-
+        battery_power_out_mw = hybrid_plant.battery.outputs.P 
+        output.update({"battery discharge [kW]": [(int(p>0))*p*1E3 for p in battery_power_out_mw]}) # convert from MW to kW and extract only discharging
+        output.update({"battery charge [kW]": [-(int(p<0))*p*1E3 for p in battery_power_out_mw]}) # convert from MW to kW and extract only charging
+        output.update({"battery state of charge [%]": hybrid_plant.battery.outputs.dispatch_SOC})
+
+    output.update({"total renewable energy production hourly [kW]": [solver_results[0]]*simulation_length})
+    output.update({"grid energy usage hourly [kW]": [solver_results[1]]*simulation_length})
+    output.update({"desal energy hourly [kW]": [solver_results[2]]*simulation_length})
+    output.update({"electrolyzer energy hourly [kW]": electrolyzer_physics_results["power_to_electrolyzer_kw"]})
+    output.update({"transport compressor energy hourly [kW]": [solver_results[3]]*simulation_length})
+    output.update({"storage energy hourly [kW]": [solver_results[4]]*simulation_length})
+    output.update({"h2 production hourly [kg]": electrolyzer_physics_results["H2_Results"]["Hydrogen Hourly Production [kg/hr]"]})
+
     df = pd.DataFrame.from_dict(output)
 
     filepath = os.path.abspath(output_dir + "data/production/")
+
     if not os.path.exists(filepath):
         os.makedirs(filepath)
 
-    df.to_csv(os.path.join(filepath, "power_series.csv"))
+    df.to_csv(os.path.join(filepath, "energy_flows.csv"))
 
-    return 0
+    return output
 
 
 # set up function to post-process HOPP results
@@ -1576,6 +1589,10 @@ def post_process_simulation(
         )
 
     # save production information
-    save_power_series(hopp_results["hybrid_plant"])
+    hourly_energy_breakdown = save_energy_flows(hopp_results["hybrid_plant"], electrolyzer_physics_results, solver_results, hours)
+
+    # save hydrogen information
+    key = "Hydrogen Hourly Production [kg/hr]"
+    np.savetxt(output_dir+"h2_usage", electrolyzer_physics_results["H2_Results"][key], header="# "+key)
 
-    return annual_energy_breakdown
+    return annual_energy_breakdown, hourly_energy_breakdown

greenheart/tools/optimization/gc_PoseOptimization.py

@@ -388,6 +388,12 @@ def set_constraints(self, opt_prob, hi: Optional[Union[None, HoppInterface]] = N
             opt_prob.model.add_subsystem("con_boundary", subsys=BoundaryDistanceComponent(hopp_interface=self.config.greenheart_config, turbine_x_init=turbine_x_init, turbine_y_init=turbine_y_init), promotes=["*"])
             opt_prob.model.add_constraint("boundary_distance_vec", lower=0)
 
+        # solar/platform size
+        # if self.config.greenheart_config["opt_options"]["constraints"]["solar_platform_ratio"]["flag"]:
+        #     upper = self.config.greenheart_config["opt_options"]["constraints"]["solar_platform_ratio"]["upper"]
+        #     opt_prob.model.add_subsystem("con_boundary", subsys=BoundaryDistanceComponent(hopp_interface=self.config.greenheart_config, turbine_x_init=turbine_x_init, turbine_y_init=turbine_y_init), promotes=["*"])
+        #     opt_prob.model.add_constraint("boundary_distance_vec", lower=0)
+
         # User constraints
         user_constr = self.config.greenheart_config["opt_options"]["constraints"]["user"]
         for k in range(len(user_constr)):

greenheart/tools/optimization/openmdao.py

@@ -26,6 +26,7 @@ def initialize(self):
 
     def setup(self):
         ninputs = 0
+        # Add inputs
         if "turbine_x" in self.options["design_variables"]:
            self.add_input("turbine_x", val=self.options["turbine_x_init"], units="m")
            ninputs += len(self.options["turbine_x_init"])
@@ -48,12 +49,23 @@ def setup(self):
             self.add_input("electrolyzer_rating_kw", val=self.options["config"].greenheart_config["electrolyzer"]["rating"]*1E3, units="kW")
             ninputs += 1
 
+        # add outputs
         self.add_output("lcoe", units="USD/(kW*h)", val=0.0, desc="levelized cost of energy")
         self.add_output("lcoh", units="USD/kg", val=0.0, desc="levelized cost of hydrogen")
         if "steel" in self.options["config"].greenheart_config.keys():
             self.add_output("lcos", units="USD/t", val=0.0, desc="levelized cost of steel")
         if "ammonia" in self.options["config"].greenheart_config.keys():
             self.add_output("lcoa", units="USD/kg", val=0.0, desc="levelized cost of ammonia")
+        if "pv_capacity_kw" in self.options["design_variables"]:
+            design_scenario = self.options["config"].plant_design_scenario
+            if self.options["config"]["greenheart_config"]["design_scenarios"][design_scenario]["pv_location"] == "offshore":
+                if self.options["config"]["greenheart_config"]["opt_options"]["constraints"]["solar_platform_ratio"]["flag"]:
+
+                    self.add_output("solar_area", units="m^2", val=0.0, desc="offshore pv array area")
+                    self.add_output("platform_area", units="m^2", val=0.0, desc="offshore platform area")
+
+                    self.options["outputs_for_finite_difference"].append(["solar_area"])
+                    self.options["outputs_for_finite_difference"].append(["platform_area"])
 
     def compute(self, inputs, outputs):
 
@@ -90,6 +102,14 @@ def compute(self, inputs, outputs):
             hopp_results, electrolyzer_physics_results, remaining_power_profile = run_simulation(config)
         elif config.output_level == 7:
             lcoe, lcoh, steel_finance, ammonia_finance = run_simulation(config)
+        elif config.output_level == 8:
+            greenheart_output = run_simulation(config)
+            lcoe = greenheart_output.lcoe
+            lcoh = greenheart_output.lcoh
+            steel_finance = greenheart_output.steel_finance
+            ammonia_finance = greenheart_output.ammonia_finance
+            # solar_area = greenheart_output.
+            # platform_area = greenheart_output
 
         outputs["lcoe"] = lcoe
         outputs["lcoh"] = lcoh