diff --git a/docs/misc_resources/turbine_models_library_preprocessing.ipynb b/docs/misc_resources/turbine_models_library_preprocessing.ipynb index 7cf4a5267..b8819ee7f 100644 --- a/docs/misc_resources/turbine_models_library_preprocessing.ipynb +++ b/docs/misc_resources/turbine_models_library_preprocessing.ipynb @@ -28,9 +28,31 @@ }, { "cell_type": "code", - "execution_count": 8, + "execution_count": 1, "metadata": {}, - "outputs": [], + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "/Users/egrant/opt/anaconda3/envs/h2i_env/lib/python3.11/site-packages/polars/_cpu_check.py:250: RuntimeWarning: Missing required CPU features.\n", + "\n", + "The following required CPU features were not detected:\n", + " avx, avx2, fma, bmi1, bmi2, lzcnt, movbe\n", + "Continuing to use this version of Polars on this processor will likely result in a crash.\n", + "Install `polars[rtcompat]` instead of `polars` to run Polars with better compatibility.\n", + "\n", + "Hint: If you are on an Apple ARM machine (e.g. M1) this is likely due to running Python under Rosetta.\n", + "It is recommended to install a native version of Python that does not run under Rosetta x86-64 emulation.\n", + "\n", + "If you believe this warning to be a false positive, you can set the `POLARS_SKIP_CPU_CHECK` environment variable to bypass this check.\n", + "\n", + " warnings.warn(\n", + "/Users/egrant/opt/anaconda3/envs/h2i_env/lib/python3.11/site-packages/fastkml/__init__.py:28: UserWarning: pkg_resources is deprecated as an API. See https://setuptools.pypa.io/en/latest/pkg_resources.html. The pkg_resources package is slated for removal as early as 2025-11-30. Refrain from using this package or pin to Setuptools<81.\n", + " from pkg_resources import DistributionNotFound\n" + ] + } + ], "source": [ "import os\n", "import numpy as np\n", @@ -52,7 +74,7 @@ }, { "cell_type": "code", - "execution_count": 9, + "execution_count": 2, "metadata": {}, "outputs": [], "source": [ @@ -104,14 +126,14 @@ }, { "cell_type": "code", - "execution_count": 10, + "execution_count": 3, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ - "/Users/gstarke/Documents/Research_Programs/H2I/H2Integrate/library/pysam_options_NREL_5MW.yaml\n" + "/Users/egrant/Documents/projects/GreenHEART/library/pysam_options_NREL_5MW.yaml\n" ] } ], @@ -125,7 +147,7 @@ }, { "cell_type": "code", - "execution_count": 11, + "execution_count": 4, "metadata": {}, "outputs": [ { @@ -367,7 +389,7 @@ " 'wind_turbine_hub_ht': 90}}" ] }, - "execution_count": 11, + "execution_count": 4, "metadata": {}, "output_type": "execute_result" } @@ -380,7 +402,7 @@ }, { "cell_type": "code", - "execution_count": 12, + "execution_count": 5, "metadata": {}, "outputs": [ { @@ -633,7 +655,7 @@ " 'layout_shape': 'square'}}}" ] }, - "execution_count": 12, + "execution_count": 5, "metadata": {}, "output_type": "execute_result" } @@ -666,7 +688,7 @@ }, { "cell_type": "code", - "execution_count": 13, + "execution_count": 6, "metadata": {}, "outputs": [ { @@ -712,7 +734,7 @@ }, { "cell_type": "code", - "execution_count": 14, + "execution_count": 7, "metadata": {}, "outputs": [ { @@ -774,7 +796,337 @@ "source": [ "## Turbine Model Pre-Processing with FLORIS\n", "\n", - "The function `export_turbine_to_floris_format()` will save turbine model specifications formatted for FLORIS. " + "Example 26 (`26_floris`) currently uses an 660 kW turbine. This example uses the \"floris_wind_plant_performance\" performance model for the wind plant. Currently, the performance model is using an 660 kW wind turbine with a rotor diameter of 47.0 meters and a hub-height of 65 meters. In the following sections we will demonstrate how to:\n", + "\n", + "1. Save turbine model specifications for the Vestas 1.65 MW turbine in the FLORIS format using `export_turbine_to_floris_format()`\n", + "2. Load the turbine model specifications for the Vestas 1.65 MW turbine and update performance parameters for the wind technology in the `tech_config` dictionary for the Vestas 1.65 MW turbine.\n", + "3. Run H2I with the updated tech_config dictionary for the Vestas 1.65 MW turbine\n", + "\n", + "\n", + "We'll start off with Step 1 and importing the function `export_turbine_to_floris_format()`, which will save turbine model specifications of the Vestas 1.65 MW turbine formatted for FLORIS. " + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "/Users/egrant/Documents/projects/GreenHEART/library/floris_turbine_Vestas_1.65MW.yaml\n" + ] + } + ], + "source": [ + "from h2integrate.preprocess.wind_turbine_file_tools import export_turbine_to_floris_format\n", + "\n", + "turbine_name = \"Vestas_1.65MW\"\n", + "\n", + "turbine_model_fpath = export_turbine_to_floris_format(turbine_name)\n", + "\n", + "print(turbine_model_fpath)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Step 2: Load the turbine model specifications for the Vestas 1.65 MW turbine and update performance parameters for the wind technology in the `tech_config` dictionary for the Vestas 1.65 MW turbine." + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "{'num_turbines': 100,\n", + " 'hub_height': -1,\n", + " 'operational_losses': 12.83,\n", + " 'floris_wake_config': {'name': 'Gauss',\n", + " 'description': 'Onshore template',\n", + " 'floris_version': 'v4.0.0',\n", + " 'logging': {'console': {'enable': False, 'level': 'WARNING'},\n", + " 'file': {'enable': False, 'level': 'WARNING'}},\n", + " 'solver': {'type': 'turbine_grid', 'turbine_grid_points': 1},\n", + " 'flow_field': {'air_density': 1.225,\n", + " 'reference_wind_height': -1,\n", + " 'wind_directions': [],\n", + " 'wind_shear': 0.33,\n", + " 'wind_speeds': [],\n", + " 'wind_veer': 0.0},\n", + " 'wake': {'model_strings': {'combination_model': 'sosfs',\n", + " 'deflection_model': 'gauss',\n", + " 'turbulence_model': 'crespo_hernandez',\n", + " 'velocity_model': 'gauss'},\n", + " 'enable_secondary_steering': False,\n", + " 'enable_yaw_added_recovery': False,\n", + " 'enable_transverse_velocities': False,\n", + " 'wake_deflection_parameters': {'gauss': {'ad': 0.0,\n", + " 'alpha': 0.58,\n", + " 'bd': 0.0,\n", + " 'beta': 0.077,\n", + " 'dm': 1.0,\n", + " 'ka': 0.38,\n", + " 'kb': 0.004},\n", + " 'jimenez': {'ad': 0.0, 'bd': 0.0, 'kd': 0.05}},\n", + " 'wake_velocity_parameters': {'cc': {'a_f': 3.11,\n", + " 'a_s': 0.179367259,\n", + " 'alpha_mod': 1.0,\n", + " 'b_f': -0.68,\n", + " 'b_s': 0.0118889215,\n", + " 'c_f': 2.41,\n", + " 'c_s1': 0.0563691592,\n", + " 'c_s2': 0.13290157},\n", + " 'gauss': {'alpha': 0.58, 'beta': 0.077, 'ka': 0.38, 'kb': 0.004},\n", + " 'jensen': {'we': 0.05}},\n", + " 'wake_turbulence_parameters': {'crespo_hernandez': {'initial': 0.1,\n", + " 'constant': 0.5,\n", + " 'ai': 0.8,\n", + " 'downstream': -0.32}},\n", + " 'enable_active_wake_mixing': False},\n", + " 'farm': {'layout_x': [], 'layout_y': []}},\n", + " 'floris_turbine_config': {'turbine_type': 'Vestas_1.65MW',\n", + " 'hub_height': 70.0,\n", + " 'TSR': 8.0,\n", + " 'rotor_diameter': 82,\n", + " 'power_thrust_table': {'wind_speed': [0.0,\n", + " 1.0,\n", + " 2.0,\n", + " 3.0,\n", + " 4.0,\n", + " 5.0,\n", + " 6.0,\n", + " 7.0,\n", + " 8.0,\n", + " 9.0,\n", + " 10.0,\n", + " 11.0,\n", + " 12.0,\n", + " 13.0,\n", + " 14.0,\n", + " 15.0,\n", + " 16.0,\n", + " 17.0,\n", + " 18.0,\n", + " 19.0,\n", + " 20.0,\n", + " 21.0,\n", + " 22.0,\n", + " 23.0,\n", + " 24.0,\n", + " 25.0,\n", + " 26.0,\n", + " 27.0,\n", + " 28.0,\n", + " 29.0,\n", + " 30.0,\n", + " 31.0,\n", + " 32.0,\n", + " 33.0,\n", + " 34.0,\n", + " 35.0,\n", + " 36.0,\n", + " 37.0,\n", + " 38.0,\n", + " 39.0,\n", + " 40.0,\n", + " 41.0,\n", + " 42.0,\n", + " 43.0,\n", + " 44.0,\n", + " 45.0,\n", + " 46.0,\n", + " 47.0,\n", + " 48.0,\n", + " 49.0],\n", + " 'power': [0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 28.0,\n", + " 144.0,\n", + " 309.0,\n", + " 511.0,\n", + " 758.0,\n", + " 1017.0,\n", + " 1285.0,\n", + " 1504.0,\n", + " 1637.0,\n", + " 1650.0,\n", + " 1650.0,\n", + " 1650.0,\n", + " 1650.0,\n", + " 1650.0,\n", + " 1650.0,\n", + " 1650.0,\n", + " 1650.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0],\n", + " 'thrust_coefficient': [0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.979,\n", + " 1.111,\n", + " 1.014,\n", + " 0.925,\n", + " 0.843,\n", + " 0.768,\n", + " 0.701,\n", + " 0.642,\n", + " 0.578,\n", + " 0.509,\n", + " 0.438,\n", + " 0.379,\n", + " 0.334,\n", + " 0.299,\n", + " 0.272,\n", + " 0.249,\n", + " 0.232,\n", + " 0.218,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0,\n", + " 0.0],\n", + " 'ref_air_density': 1.225,\n", + " 'ref_tilt': 5.0,\n", + " 'cosine_loss_exponent_yaw': 1.88,\n", + " 'cosine_loss_exponent_tilt': 1.88}},\n", + " 'resource_data_averaging_method': 'weighted_average',\n", + " 'floris_operation_model': 'cosine-loss',\n", + " 'default_turbulence_intensity': 0.06,\n", + " 'enable_caching': True,\n", + " 'cache_dir': 'cache',\n", + " 'layout': {'layout_mode': 'basicgrid',\n", + " 'layout_options': {'row_D_spacing': 5.0,\n", + " 'turbine_D_spacing': 5.0,\n", + " 'rotation_angle_deg': 0.0,\n", + " 'row_phase_offset': 0.0,\n", + " 'layout_shape': 'square'}}}" + ] + }, + "execution_count": 9, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "# Load the tech config file\n", + "tech_config_path = EXAMPLE_DIR / \"26_floris\" / \"tech_config.yaml\"\n", + "tech_config = load_tech_yaml(tech_config_path)\n", + "\n", + "# Load the turbine model file formatted for FLORIS\n", + "floris_options = load_yaml(turbine_model_fpath)\n", + "\n", + "# Create dictionary of updated inputs for the new turbine formatted for\n", + "# the \"floris_wind_plant_performance\" model\n", + "updated_parameters = {\n", + " \"hub_height\": -1, # -1 indicates to use the hub-height in the floris_turbine_config\n", + " \"floris_turbine_config\": floris_options,\n", + "}\n", + "\n", + "# Update wind performance parameters in the tech config\n", + "tech_config[\"technologies\"][\"wind\"][\"model_inputs\"][\"performance_parameters\"].update(\n", + " updated_parameters\n", + ")\n", + "\n", + "# The technology input for the updated wind turbine model\n", + "tech_config[\"technologies\"][\"wind\"][\"model_inputs\"][\"performance_parameters\"]" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Step 3: Run H2I with the updated tech_config dictionary for the Vestas 1.65 MW turbine" + ] + }, + { + "cell_type": "code", + "execution_count": 11, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Wind LCOE is $120.87/MWh\n" + ] + } + ], + "source": [ + "# Create the top-level config input dictionary for H2I\n", + "h2i_config = {\n", + " \"driver_config\": EXAMPLE_DIR / \"26_floris\" / \"driver_config.yaml\",\n", + " \"technology_config\": tech_config,\n", + " \"plant_config\": EXAMPLE_DIR / \"26_floris\" / \"plant_config.yaml\",\n", + "}\n", + "\n", + "# Create a H2Integrate model with the updated tech config\n", + "h2i = H2IntegrateModel(h2i_config)\n", + "\n", + "# Run the model\n", + "h2i.run()\n", + "\n", + "# Get LCOE of wind plant\n", + "wind_lcoe = h2i.model.get_val(\"finance_subgroup_electricity.LCOE\", units=\"USD/MW/h\")\n", + "print(f\"Wind LCOE is ${wind_lcoe[0]:.2f}/MWh\")" ] } ], diff --git a/docs/user_guide/model_overview.md b/docs/user_guide/model_overview.md index f79c277ab..03adfa715 100644 --- a/docs/user_guide/model_overview.md +++ b/docs/user_guide/model_overview.md @@ -133,6 +133,7 @@ Below summarizes the available performance, cost, and financial models for each - `wind`: wind turbine - performance models: + `'pysam_wind_plant_performance'` + + `'floris_wind_plant_performance'` - cost models: + `'atb_wind_cost'` - `solar`: solar-PV panels diff --git a/examples/26_floris/driver_config.yaml b/examples/26_floris/driver_config.yaml new file mode 100644 index 000000000..aff911d3b --- /dev/null +++ b/examples/26_floris/driver_config.yaml @@ -0,0 +1,5 @@ +name: "driver_config" +description: "This analysis runs a wind plant using FLORIS" + +general: + folder_output: floris_wind_outputs diff --git a/examples/26_floris/plant_config.yaml b/examples/26_floris/plant_config.yaml new file mode 100644 index 000000000..a12072147 --- /dev/null +++ b/examples/26_floris/plant_config.yaml @@ -0,0 +1,54 @@ +name: "plant_config" +description: "" + +site: + latitude: 44.04218 + longitude: -95.19757 + + resources: + wind_resource: + resource_model: "openmeteo_wind_api" + resource_parameters: + resource_year: 2023 + +resource_to_tech_connections: [ + # connect the wind resource to the wind technology + ['wind_resource', 'wind', 'wind_resource_data'], +] + +plant: + plant_life: 30 + simulation: + n_timesteps: 8760 + dt: 3600 +finance_parameters: + finance_groups: + finance_model: "ProFastComp" + model_inputs: + params: + analysis_start_year: 2032 + installation_time: 36 # months + inflation_rate: 0.0 # 0 for nominal analysis + discount_rate: 0.09 + debt_equity_ratio: 2.62 + property_tax_and_insurance: 0.03 # percent of CAPEX + total_income_tax_rate: 0.257 + capital_gains_tax_rate: 0.15 # H2FAST default + sales_tax_rate: 0.07375 + debt_interest_rate: 0.07 + debt_type: "Revolving debt" # can be "Revolving debt" or "One time loan". + 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 #years + refurb: [0.] + finance_subgroups: + electricity: + commodity: "electricity" + commodity_stream: "wind" #use the total electricity output from the combiner for the finance calc + technologies: ["wind"] + cost_adjustment_parameters: + cost_year_adjustment_inflation: 0.025 # used to adjust modeled costs to target_dollar_year + target_dollar_year: 2022 diff --git a/examples/26_floris/run_floris_example.py b/examples/26_floris/run_floris_example.py new file mode 100644 index 000000000..1a39a494b --- /dev/null +++ b/examples/26_floris/run_floris_example.py @@ -0,0 +1,27 @@ +import os +from pathlib import Path + +from h2integrate.core.utilities import load_yaml +from h2integrate.core.h2integrate_model import H2IntegrateModel + + +os.chdir(Path(__file__).parent) + +driver_config = load_yaml(Path(__file__).parent / "driver_config.yaml") +tech_config = load_yaml(Path(__file__).parent / "tech_config.yaml") +plant_config = load_yaml(Path(__file__).parent / "plant_config.yaml") + +h2i_config = { + "name": "H2Integrate_config", + "system_summary": "", + "driver_config": driver_config, + "technology_config": tech_config, + "plant_config": plant_config, +} +# Create a H2I model +h2i = H2IntegrateModel(h2i_config) + +# Run the model +h2i.run() + +h2i.post_process() diff --git a/examples/26_floris/tech_config.yaml b/examples/26_floris/tech_config.yaml new file mode 100644 index 000000000..833bf636d --- /dev/null +++ b/examples/26_floris/tech_config.yaml @@ -0,0 +1,33 @@ +name: "technology_config" +description: "This plant runs a wind farm using FLORIS" + +technologies: + wind: + performance_model: + model: "floris_wind_plant_performance" + cost_model: + model: "atb_wind_cost" + model_inputs: + performance_parameters: + num_turbines: 100 #number of turbines in the farm + hub_height: 65.0 # turbine hub-height + operational_losses: 12.83 #percentage of non-wake losses + floris_wake_config: !include "floris_v4_default_template.yaml" #floris wake model file + floris_turbine_config: !include "floris_turbine_Vestas_660kW.yaml" #turbine model file formatted for floris + resource_data_averaging_method: "weighted_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: 5.0 + turbine_D_spacing: 5.0 + rotation_angle_deg: 0.0 + row_phase_offset: 0.0 + layout_shape: "square" + cost_parameters: + capex_per_kW: 1472.0 + opex_per_kW_per_year: 32 + cost_year: 2022 diff --git a/examples/test/test_all_examples.py b/examples/test/test_all_examples.py index d8de2fbe3..18e51e5c4 100644 --- a/examples/test/test_all_examples.py +++ b/examples/test/test_all_examples.py @@ -1260,6 +1260,56 @@ def test_sweeping_solar_sites_doe(subtests): assert len(list(set(res_df["LCOE"].to_list()))) == len(res_df) +def test_floris_example(subtests): + from h2integrate.core.utilities import load_yaml + + os.chdir(EXAMPLE_DIR / "26_floris") + + driver_config = load_yaml(EXAMPLE_DIR / "26_floris" / "driver_config.yaml") + tech_config = load_yaml(EXAMPLE_DIR / "26_floris" / "tech_config.yaml") + plant_config = load_yaml(EXAMPLE_DIR / "26_floris" / "plant_config.yaml") + + h2i_config = { + "name": "H2Integrate_config", + "system_summary": "", + "driver_config": driver_config, + "technology_config": tech_config, + "plant_config": plant_config, + } + + # Create a H2I model + h2i = H2IntegrateModel(h2i_config) + + # Run the model + h2i.run() + + with subtests.test("LCOE"): + assert ( + pytest.approx( + h2i.prob.get_val("finance_subgroup_electricity.LCOE", units="USD/MW/h")[0], rel=1e-6 + ) + == 99.872209 + ) + + with subtests.test("Wind plant capacity"): + assert pytest.approx(h2i.prob.get_val("wind.total_capacity", units="MW"), rel=1e-6) == 66.0 + + with subtests.test("Total electricity production"): + assert ( + pytest.approx( + np.sum(h2i.prob.get_val("wind.total_electricity_produced", units="MW*h/yr")), + rel=1e-6, + ) + == 128948.21977 + ) + + with subtests.test("Capacity factor"): + assert ( + pytest.approx(h2i.prob.get_val("wind.capacity_factor", units="percent")[0], rel=1e-6) + == 22.30320668 + ) + + def test_24_solar_battery_grid_example(subtests): # NOTE: would be good to compare LCOE against the same example without grid selling # and see that LCOE reduces with grid selling diff --git a/h2integrate/converters/hopp/test/test_hopp_caching.py b/h2integrate/converters/hopp/test/test_hopp_caching.py index 5db630aa9..2579826ed 100644 --- a/h2integrate/converters/hopp/test/test_hopp_caching.py +++ b/h2integrate/converters/hopp/test/test_hopp_caching.py @@ -72,3 +72,6 @@ def test_hopp_wrapper_cache_filenames(subtests, plant_config, tech_config): with subtests.test("Check unique filename with modified config"): assert len(cache_filename_new) > 0 + + # Delete cache files and the testing cache dir + shutil.rmtree(cache_dir) diff --git a/h2integrate/converters/wind/floris.py b/h2integrate/converters/wind/floris.py new file mode 100644 index 000000000..623b68322 --- /dev/null +++ b/h2integrate/converters/wind/floris.py @@ -0,0 +1,294 @@ +import copy + +import numpy as np +from attrs import field, define, validators +from floris import TimeSeries, FlorisModel + +from h2integrate.core.utilities import merge_shared_inputs +from h2integrate.core.validators import gt_zero, contains, range_val +from h2integrate.core.model_baseclasses import CacheBaseClass, CacheBaseConfig +from h2integrate.converters.wind.tools.resource_tools import ( + calculate_air_density, + average_wind_data_for_hubheight, + weighted_average_wind_data_for_hubheight, +) +from h2integrate.converters.wind.wind_plant_baseclass import WindPerformanceBaseClass +from h2integrate.converters.wind.layout.simple_grid_layout import ( + BasicGridLayoutConfig, + make_basic_grid_turbine_layout, +) + + +@define +class FlorisWindPlantPerformanceConfig(CacheBaseConfig): + """Configuration class for FlorisWindPlantPerformanceModel. + + Attributes: + num_turbines (int): number of turbines in farm + floris_wake_config (dict): dictionary containing FLORIS inputs for flow_field, wake, + solver, and logging. + floris_turbine_config (dict): dictionary containing turbine parameters formatted for FLORIS. + operational_losses (float | int): non-wake losses represented as a percentage + (between 0 and 100). + hub_height (float | int, optional): a value of -1 indicates to use the hub-height + from the `floris_turbine_config`. Otherwise, is the turbine hub-height + in meters. Defaults to -1. + operation_model (str, optional): turbine operation model. Defaults to 'cosine-loss'. + default_turbulence_intensity (float): default turbulence intensity to use if not found + in wind resource data. + layout (dict): layout parameters dictionary. + resource_data_averaging_method (str): string indicating what method to use to + adjust or select resource data if no resource data is available at a height + exactly equal to the turbine hub-height. Defaults to 'weighted_average'. + The available methods are: + - 'weighted_average': average the resource data at the heights that most closely bound + the hub-height, weighted by the difference between the resource heights and the + hub-height. + - 'average': average the resource data at the heights that most closely bound + the hub-height. + - 'nearest': use the resource data at the height closest to the hub-height. + enable_caching (bool, optional): if True, checks if the outputs have be saved to a + cached file or saves outputs to a file. Defaults to True. + cache_dir (str | Path, optional): folder to use for reading or writing cached results files. + Only used if enable_caching is True. Defaults to "cache". + hybrid_turbine_design (bool, optional): whether multiple turbine types are included in + the farm. Defaults to False. The functionality to use multiple turbine types + is not yet implemented. Will result in NotImplementedError if True. + """ + + num_turbines: int = field(converter=int, validator=gt_zero) + floris_wake_config: dict = field() + floris_turbine_config: dict = field() + default_turbulence_intensity: float = field() + operational_losses: float = field(validator=range_val(0.0, 100.0)) + hub_height: float = field(default=-1, validator=validators.ge(-1)) + adjust_air_density_for_elevation: bool = field(default=False) + operation_model: str = field(default="cosine-loss") + layout: dict = field(default={}) + resource_data_averaging_method: str = field( + default="weighted_average", validator=contains(["weighted_average", "average", "nearest"]) + ) + hybrid_turbine_design: bool = field(default=False) + + # if using multiple turbines, then need to specify resource reference height + def __attrs_post_init__(self): + super().__attrs_post_init__() + n_turbine_types = len(self.floris_wake_config.get("farm", {}).get("turbine_type", [])) + n_pos = len(self.floris_wake_config.get("farm", {}).get("layout_x", [])) + if n_turbine_types > 1 and n_turbine_types != n_pos: + self.hybrid_turbine_design = True + + # use floris_turbines + if self.hub_height < 0 and not self.hybrid_turbine_design: + self.hub_height = self.floris_turbine_config.get("hub_height") + + # check that user did not provide a layout in the floris_wake_config + gave_x_coords = len(self.floris_wake_config.get("farm", {}).get("layout_x", [])) > 0 + gave_y_coords = len(self.floris_wake_config.get("farm", {}).get("layout_y", [])) > 0 + if gave_x_coords or gave_y_coords: + msg = ( + "Layout provided in `floris_wake_config['farm']` but layout will be created " + "based on the `layout_mode` and `layout_options` provided in the " + "`layout` dictionary. Please set the layout in " + "floris_wake_config['farm']['layout_x'] and floris_wake_config['farm']['layout_y']" + " to empty lists" + ) + raise ValueError(msg) + + +class FlorisWindPlantPerformanceModel(WindPerformanceBaseClass, CacheBaseClass): + """ + An OpenMDAO component that wraps a Floris model. + It takes wind turbine model parameters and wind resource data as input and + outputs power generation data. + """ + + def setup(self): + self.n_timesteps = int(self.options["plant_config"]["plant"]["simulation"]["n_timesteps"]) + + performance_inputs = self.options["tech_config"]["model_inputs"]["performance_parameters"] + + # initialize layout config + layout_options = {} + if "layout" in performance_inputs: + layout_params = self.options["tech_config"]["model_inputs"]["performance_parameters"][ + "layout" + ] + layout_mode = layout_params.get("layout_mode", "basicgrid") + layout_options = layout_params.get("layout_options", {}) + if layout_mode == "basicgrid": + self.layout_config = BasicGridLayoutConfig.from_dict(layout_options) + self.layout_mode = layout_mode + + # initialize wind turbine config + self.config = FlorisWindPlantPerformanceConfig.from_dict( + merge_shared_inputs(self.options["tech_config"]["model_inputs"], "performance") + ) + + if self.config.hybrid_turbine_design: + raise NotImplementedError( + "H2I does not currently support running multiple different wind turbine " + "designs with Floris." + ) + + self.add_input( + "num_turbines", + val=self.config.num_turbines, + units="unitless", + desc="number of turbines in farm", + ) + + self.add_input( + "hub_height", + val=self.config.hub_height, + units="m", + desc="turbine hub-height", + ) + + self.add_output( + "total_electricity_produced", + val=0.0, + units="kW*h/year", + desc="Annual energy production from WindPlant", + ) + self.add_output("total_capacity", val=0.0, units="kW", desc="Wind farm rated capacity") + + self.add_output( + "capacity_factor", val=0.0, units="unitless", desc="Wind farm capacity factor" + ) + + super().setup() + + power_curve = self.config.floris_turbine_config.get("power_thrust_table").get("power") + self.wind_turbine_rating_kW = np.max(power_curve) + + def format_resource_data(self, hub_height, wind_resource_data): + # NOTE: could weight resource data of bounding heights like + # `weighted_parse_resource_data` method in HOPP + + bounding_heights = self.calculate_bounding_heights_from_resource_data( + hub_height, wind_resource_data, resource_vars=["wind_speed", "wind_direction"] + ) + if len(bounding_heights) == 1: + resource_height = bounding_heights[0] + windspeed = wind_resource_data[f"wind_speed_{resource_height}m"] + winddir = wind_resource_data[f"wind_direction_{resource_height}m"] + else: + if self.config.resource_data_averaging_method == "nearest": + height_difference = [np.abs(hub_height - b) for b in bounding_heights] + resource_height = bounding_heights[np.argmin(height_difference).flatten()[0]] + windspeed = wind_resource_data[f"wind_speed_{resource_height}m"] + winddir = wind_resource_data[f"wind_direction_{resource_height}m"] + if self.config.resource_data_averaging_method == "weighted_average": + windspeed = weighted_average_wind_data_for_hubheight( + wind_resource_data, bounding_heights, hub_height, "wind_speed" + ) + winddir = weighted_average_wind_data_for_hubheight( + wind_resource_data, bounding_heights, hub_height, "wind_direction" + ) + if self.config.resource_data_averaging_method == "average": + windspeed = average_wind_data_for_hubheight( + wind_resource_data, bounding_heights, "wind_speed" + ) + winddir = average_wind_data_for_hubheight( + wind_resource_data, bounding_heights, "wind_direction" + ) + + # get turbulence intensity + + # check if turbulence intensity is available in wind resource data + if any("turbulence_intensity_" in k for k in wind_resource_data.keys()): + for height in bounding_heights: + if f"turbulence_intensity_{height}m" in wind_resource_data: + ti = wind_resource_data.get( + f"turbulence_intensity_{height}m", self.config.default_turbulence_intensity + ) + break + else: + ti = wind_resource_data.get( + "turbulence_intensity", self.config.default_turbulence_intensity + ) + + time_series = TimeSeries( + wind_directions=winddir, + wind_speeds=windspeed, + turbulence_intensities=ti, + ) + + return time_series + + def compute(self, inputs, outputs, discrete_inputs, discrete_outputs): + # 1. Check if the results for the current configuration are already cached + config_dict = self.config.as_dict() + config_dict.update({"wind_turbine_size_kw": self.wind_turbine_rating_kW}) + loaded_results = self.load_outputs( + inputs, outputs, discrete_inputs, discrete_outputs={}, config_dict=config_dict + ) + if loaded_results: + # Case has been run before and outputs have been set, can exit this function + return + + # 2. If caching is not enabled or a cache file does not exist, run FLORIS + n_turbs = int(np.round(inputs["num_turbines"][0])) + + # Copy main config files + floris_config = copy.deepcopy(self.config.floris_wake_config) + turbine_design = copy.deepcopy(self.config.floris_turbine_config) + + # update the turbine hub-height in the floris turbine config + turbine_design.update({"hub_height": inputs["hub_height"][0]}) + + # update the operation model in the floris turbine config + turbine_design.update({"operation_model": self.config.operation_model}) + + # format resource data and input into model + time_series = self.format_resource_data( + inputs["hub_height"][0], discrete_inputs["wind_resource_data"] + ) + + # make layout for number of turbines + if self.layout_mode == "basicgrid": + x_pos, y_pos = make_basic_grid_turbine_layout( + turbine_design.get("rotor_diameter"), n_turbs, self.layout_config + ) + + floris_farm = {"layout_x": x_pos, "layout_y": y_pos, "turbine_type": [turbine_design]} + + floris_config["farm"].update(floris_farm) + + # adjust air density + if ( + self.config.adjust_air_density_for_elevation + and "elevation" in discrete_inputs["wind_resource_data"] + ): + rho = calculate_air_density(discrete_inputs["wind_resource_data"]["elevation"]) + floris_config["flow_field"].update({"air_density": rho}) + + # initialize FLORIS + floris_config["flow_field"].update({"turbulence_intensities": []}) + self.fi = FlorisModel(floris_config) + + # set the layout and wind data in Floris + self.fi.set(layout_x=x_pos, layout_y=y_pos, wind_data=time_series) + + # run the model + self.fi.run() + + power_farm = self.fi.get_farm_power().reshape(self.n_timesteps) # W + + # Adding losses (excluding turbine and wake losses) + operational_efficiency = (100 - self.config.operational_losses) / 100 + gen = power_farm * operational_efficiency / 1000 # kW + + # set outputs + outputs["electricity_out"] = gen + outputs["total_capacity"] = n_turbs * self.wind_turbine_rating_kW + + max_production = n_turbs * self.wind_turbine_rating_kW * len(gen) + outputs["total_electricity_produced"] = np.sum(gen) + outputs["capacity_factor"] = np.sum(gen) / max_production + + # 3. Cache the results for future use if enabled + self.cache_outputs( + inputs, outputs, discrete_inputs, discrete_outputs={}, config_dict=config_dict + ) diff --git a/h2integrate/converters/wind/test/test_floris_wind.py b/h2integrate/converters/wind/test/test_floris_wind.py new file mode 100644 index 000000000..4a0e0d0d5 --- /dev/null +++ b/h2integrate/converters/wind/test/test_floris_wind.py @@ -0,0 +1,364 @@ +import shutil + +import numpy as np +import pytest +import openmdao.api as om +from pytest import fixture + +from h2integrate import ROOT_DIR, H2I_LIBRARY_DIR +from h2integrate.core.utilities import load_yaml +from h2integrate.converters.wind.floris import FlorisWindPlantPerformanceModel +from h2integrate.resource.wind.openmeteo_wind import OpenMeteoHistoricalWindResource +from h2integrate.resource.wind.nrel_developer_wtk_api import WTKNRELDeveloperAPIWindResource + + +@fixture +def floris_config(): + floris_wake_config = load_yaml(H2I_LIBRARY_DIR / "floris_v4_default_template.yaml") + floris_turbine_config = load_yaml(H2I_LIBRARY_DIR / "floris_turbine_Vestas_660kW.yaml") + floris_performance_dict = { + "num_turbines": 20, + "floris_wake_config": floris_wake_config, + "floris_turbine_config": floris_turbine_config, + "default_turbulence_intensity": 0.06, + "operation_model": "cosine-loss", + "hub_height": -1, + "layout": { + "layout_mode": "basicgrid", + "layout_options": { + "row_D_spacing": 5.0, + "turbine_D_spacing": 5.0, + }, + }, + "operational_losses": 12.83, + "enable_caching": False, + "cache_dir": None, + "resource_data_averaging_method": "nearest", + } + return floris_performance_dict + + +@fixture +def plant_config_openmeteo(): + site_config = { + "latitude": 44.04218, + "longitude": -95.19757, + "resource": { + "wind_resource": { + "resource_model": "openmeteo_wind_api", + "resource_parameters": { + "resource_year": 2023, + }, + } + }, + } + plant_dict = { + "plant_life": 30, + "simulation": {"n_timesteps": 8760, "dt": 3600, "start_time": "01/01 00:30:00"}, + } + + d = {"site": site_config, "plant": plant_dict} + return d + + +@fixture +def plant_config_wtk(): + site_config = { + "latitude": 35.2018863, + "longitude": -101.945027, + "resource": { + "wind_resource": { + "resource_model": "wind_toolkit_v2_api", + "resource_parameters": { + "resource_year": 2012, + }, + } + }, + } + plant_dict = { + "plant_life": 30, + "simulation": {"n_timesteps": 8760, "dt": 3600, "start_time": "01/01 00:30:00"}, + } + + d = {"site": site_config, "plant": plant_dict} + return d + + +def test_floris_wind_performance(plant_config_openmeteo, floris_config, subtests): + tech_config_dict = { + "model_inputs": { + "performance_parameters": floris_config, + } + } + + prob = om.Problem() + + wind_resource_config = plant_config_openmeteo["site"]["resource"]["wind_resource"][ + "resource_parameters" + ] + wind_resource = OpenMeteoHistoricalWindResource( + plant_config=plant_config_openmeteo, + resource_config=wind_resource_config, + driver_config={}, + ) + + wind_plant = FlorisWindPlantPerformanceModel( + plant_config=plant_config_openmeteo, + tech_config=tech_config_dict, + driver_config={}, + ) + + prob.model.add_subsystem("wind_resource", wind_resource, promotes=["*"]) + prob.model.add_subsystem("wind_plant", wind_plant, promotes=["*"]) + prob.setup() + prob.run_model() + + with subtests.test("wind farm capacity"): + assert ( + pytest.approx(prob.get_val("wind_plant.total_capacity", units="kW")[0], rel=1e-6) + == 660 * 20 + ) + + with subtests.test("AEP"): + assert ( + pytest.approx( + prob.get_val("wind_plant.total_electricity_produced", units="kW*h/year")[0], + rel=1e-6, + ) + == 36471.03023616864 * 1e3 + ) + + with subtests.test("total electricity_out"): + assert pytest.approx( + np.sum(prob.get_val("wind_plant.electricity_out", units="kW")), rel=1e-6 + ) == prob.get_val("wind_plant.total_electricity_produced", units="kW*h/year") + + +def test_floris_caching_changed_config(plant_config_openmeteo, floris_config, subtests): + cache_dir = ROOT_DIR.parent / "test_cache_floris" + + # delete cache dir if it exists + if cache_dir.exists(): + shutil.rmtree(cache_dir) + + floris_config["enable_caching"] = True + floris_config["cache_dir"] = cache_dir + + tech_config_dict = { + "model_inputs": { + "performance_parameters": floris_config, + } + } + + # Run FLORIS and get cache filename + prob = om.Problem() + + wind_resource_config = plant_config_openmeteo["site"]["resource"]["wind_resource"][ + "resource_parameters" + ] + wind_resource = OpenMeteoHistoricalWindResource( + plant_config=plant_config_openmeteo, + resource_config=wind_resource_config, + driver_config={}, + ) + + wind_plant = FlorisWindPlantPerformanceModel( + plant_config=plant_config_openmeteo, + tech_config=tech_config_dict, + driver_config={}, + ) + + prob.model.add_subsystem("wind_resource", wind_resource, promotes=["*"]) + prob.model.add_subsystem("wind_plant", wind_plant, promotes=["*"]) + prob.setup() + prob.run_model() + + cache_filename_init = list(cache_dir.glob("*.pkl")) + + with subtests.test("Check that cache file was created"): + assert len(cache_filename_init) == 1 + + # Modify something in the config and check that cache filename is different + floris_config["operational_losses"] = 10.0 + prob = om.Problem() + wind_resource = OpenMeteoHistoricalWindResource( + plant_config=plant_config_openmeteo, + resource_config=wind_resource_config, + driver_config={}, + ) + + wind_plant = FlorisWindPlantPerformanceModel( + plant_config=plant_config_openmeteo, + tech_config=tech_config_dict, + driver_config={}, + ) + + prob.model.add_subsystem("wind_resource", wind_resource, promotes=["*"]) + prob.model.add_subsystem("wind_plant", wind_plant, promotes=["*"]) + prob.setup() + prob.run_model() + + cache_filenames = list(cache_dir.glob("*.pkl")) + cache_filename_new = [file for file in cache_filenames if file not in cache_filename_init] + + with subtests.test("Check unique filename with modified config"): + assert len(cache_filename_new) > 0 + + with subtests.test("Check two cache files exist"): + assert len(cache_filenames) == 2 + + # Delete cache files and the testing cache dir + shutil.rmtree(cache_dir) + + +def test_floris_caching_changed_inputs(plant_config_openmeteo, floris_config, subtests): + cache_dir = ROOT_DIR.parent / "test_cache_floris" + + # delete cache dir if it exists + if cache_dir.exists(): + shutil.rmtree(cache_dir) + + floris_config["enable_caching"] = True + floris_config["cache_dir"] = cache_dir + + tech_config_dict = { + "model_inputs": { + "performance_parameters": floris_config, + } + } + + # Run FLORIS and get cache filename + prob = om.Problem() + + wind_resource_config = plant_config_openmeteo["site"]["resource"]["wind_resource"][ + "resource_parameters" + ] + wind_resource = OpenMeteoHistoricalWindResource( + plant_config=plant_config_openmeteo, + resource_config=wind_resource_config, + driver_config={}, + ) + + wind_plant = FlorisWindPlantPerformanceModel( + plant_config=plant_config_openmeteo, + tech_config=tech_config_dict, + driver_config={}, + ) + + prob.model.add_subsystem("wind_resource", wind_resource, promotes=["*"]) + prob.model.add_subsystem("wind_plant", wind_plant, promotes=["*"]) + prob.setup() + prob.run_model() + + wind_resource_data = dict(prob.get_val("wind_resource.wind_resource_data")) + + cache_filename_init = list(cache_dir.glob("*.pkl")) + + with subtests.test("Check that cache file was created"): + assert len(cache_filename_init) == 1 + + # Modify the wind resource data, rerun floris, and check that a new file was created + # wind_resource_data['wind_speed_100m'][10] += 1 #this wont trigger a new cache file + wind_resource_data["site_lat"] = 44.04218107666016 + + prob = om.Problem() + + wind_plant = FlorisWindPlantPerformanceModel( + plant_config=plant_config_openmeteo, + tech_config=tech_config_dict, + driver_config={}, + ) + + prob.model.add_subsystem("wind_plant", wind_plant) + prob.setup() + prob.set_val("wind_plant.wind_resource_data", wind_resource_data) + prob.run_model() + + cache_filenames = list(cache_dir.glob("*.pkl")) + cache_filename_new = [file for file in cache_filenames if file not in cache_filename_init] + + with subtests.test("Check unique filename with modified config"): + assert len(cache_filename_new) > 0 + + with subtests.test("Check two cache files exist"): + assert len(cache_filenames) == 2 + + # Delete cache files and the testing cache dir + shutil.rmtree(cache_dir) + + +def test_floris_wind_performance_air_dens(plant_config_wtk, floris_config, subtests): + tech_config_dict = { + "model_inputs": { + "performance_parameters": floris_config, + } + } + + prob = om.Problem() + + wind_resource_config = plant_config_wtk["site"]["resource"]["wind_resource"][ + "resource_parameters" + ] + wind_resource = WTKNRELDeveloperAPIWindResource( + plant_config=plant_config_wtk, + resource_config=wind_resource_config, + driver_config={}, + ) + + wind_plant = FlorisWindPlantPerformanceModel( + plant_config=plant_config_wtk, + tech_config=tech_config_dict, + driver_config={}, + ) + + prob.model.add_subsystem("wind_resource", wind_resource, promotes=["*"]) + prob.model.add_subsystem("wind_plant", wind_plant, promotes=["*"]) + prob.setup() + prob.run_model() + + wind_resource_data = dict(prob.get_val("wind_resource.wind_resource_data")) + + initial_aep = prob.get_val("wind_plant.total_electricity_produced", units="kW*h/year")[0] + with subtests.test("wind farm capacity"): + assert ( + pytest.approx(prob.get_val("wind_plant.total_capacity", units="kW")[0], rel=1e-6) + == 660 * 20 + ) + + with subtests.test("AEP"): + assert ( + pytest.approx( + prob.get_val("wind_plant.total_electricity_produced", units="kW*h/year")[0], + rel=1e-6, + ) + == 37007.33639643173 * 1e3 + ) + + with subtests.test("total electricity_out"): + assert pytest.approx( + np.sum(prob.get_val("wind_plant.electricity_out", units="kW")), rel=1e-6 + ) == prob.get_val("wind_plant.total_electricity_produced", units="kW*h/year") + + # Add elevation to the resource data and rerun floris + floris_config["adjust_air_density_for_elevation"] = True + wind_resource_data["elevation"] = 1133.0 + + prob = om.Problem() + + wind_plant = FlorisWindPlantPerformanceModel( + plant_config=plant_config_wtk, + tech_config=tech_config_dict, + driver_config={}, + ) + + prob.model.add_subsystem("wind_plant", wind_plant) + prob.setup() + prob.set_val("wind_plant.wind_resource_data", wind_resource_data) + prob.run_model() + + adjusted_aep = prob.get_val("wind_plant.total_electricity_produced", units="kW*h/year")[0] + with subtests.test("reduced AEP with air density adjustment"): + assert adjusted_aep < initial_aep + + with subtests.test("AEP with air density adjustment"): + assert pytest.approx(adjusted_aep, rel=1e-6) == 34392.58173437373 * 1e3 diff --git a/h2integrate/converters/wind/tools/__init__.py b/h2integrate/converters/wind/tools/__init__.py new file mode 100644 index 000000000..e69de29bb diff --git a/h2integrate/converters/wind/tools/resource_tools.py b/h2integrate/converters/wind/tools/resource_tools.py new file mode 100644 index 000000000..2272812c9 --- /dev/null +++ b/h2integrate/converters/wind/tools/resource_tools.py @@ -0,0 +1,131 @@ +import numpy as np +from scipy.constants import R, g, convert_temperature + + +RHO_0 = 1.225 # Air density at sea level (kg/m3) +T_REF = 20 # Standard air temperature (Celsius) +MOLAR_MASS_AIR = 28.96 # Molar mass of air (g/mol) +LAPSE_RATE = 0.0065 # Temperature lapse rate (K/m) for 0-11000m above sea level + + +def calculate_air_density(elevation_m: float) -> float: + """ + Calculate air density based on site elevation using the Barometric formula. + + This function is based on Equation 1 from: https://en.wikipedia.org/wiki/Barometric_formula#Density_equations + Imported constants are: + + - g: acceleration due to gravity (m/s2) + - R: universal gas constant (J/mol-K) + + Args: + elevation_m (float): Elevation of site in meters + + Returns: + float: Air density in kg/m^3 at elevation of site + """ + + # Reference elevation at sea level (m) + elevation_sea_level = 0.0 + + # Convert temperature to Kelvin + T_ref_K = convert_temperature([T_REF], "C", "K")[0] + + # Exponent value used in equation below + e = g * (MOLAR_MASS_AIR / 1e3) / (R * LAPSE_RATE) + + # Calculate air density at site elevation + rho = RHO_0 * ((T_ref_K - ((elevation_m - elevation_sea_level) * LAPSE_RATE)) / T_ref_K) ** ( + e - 1 + ) + return rho + + +def weighted_average_wind_data_for_hubheight( + wind_resource_data: dict, + bounding_resource_heights: tuple[int] | list[int], + hub_height: float | int, + wind_resource_spec: str, +): + """Compute the weighted average of wind resource data at two resource heights. + + Args: + wind_resource_data (dict): dictionary of wind resource data + bounding_resource_heights (tuple[int] | list[int]): resource heights that bound the + hub-height, formatted as [lower_resource_height, upper_resource_height] + hub_height (float | int): wind turbine hub-height in meters. + wind_resource_spec (str): wind resource data key that is unique for + each hub-height. Such as `'wind_speed'` or `'wind_direction'` + + Raises: + ValueError: if f'{wind_resource_spec}_{lower_resource_height}m' or + f'{wind_resource_spec}_{upper_resource_height}m' are not found in `wind_resource_data` + + Returns: + np.ndarray: wind resource data averaged between the two bounding heights. + """ + height_lower, height_upper = bounding_resource_heights + + has_lowerbound = f"{wind_resource_spec}_{height_lower}m" in wind_resource_data + has_upperbound = f"{wind_resource_spec}_{height_upper}m" in wind_resource_data + if not has_lowerbound or not has_upperbound: + msg = ( + f"Wind resource data for {wind_resource_spec} is missing either " + f"{height_lower}m or {height_upper}m" + ) + raise ValueError(msg) + + # weight1 is the weight applied to the lower-bound height + weight1 = np.abs(height_upper - hub_height) + # weight2 is the weight applied to the upper-bound height + weight2 = np.abs(height_lower - hub_height) + + weighted_wind_resource = ( + (weight1 * wind_resource_data[f"{wind_resource_spec}_{height_lower}m"]) + + (weight2 * wind_resource_data[f"{wind_resource_spec}_{height_upper}m"]) + ) / (weight1 + weight2) + + return weighted_wind_resource + + +def average_wind_data_for_hubheight( + wind_resource_data: dict, + bounding_resource_heights: tuple[int] | list[int], + wind_resource_spec: str, +): + """Compute the average of wind resource data at two resource heights. + + Args: + wind_resource_data (dict): dictionary of wind resource data + bounding_resource_heights (tuple[int] | list[int]): resource heights that bound the + hub-height, formatted as [lower_resource_height, upper_resource_height] + wind_resource_spec (str): wind resource data key that is unique for + each hub-height. Such as `'wind_speed'` or `'wind_direction'` + + Raises: + ValueError: if f'{wind_resource_spec}_{lower_resource_height}m' or + f'{wind_resource_spec}_{upper_resource_height}m' are not found in `wind_resource_data` + + Returns: + np.ndarray: wind resource data averaged between the two bounding heights. + """ + height_lower, height_upper = bounding_resource_heights + + has_lowerbound = f"{wind_resource_spec}_{height_lower}m" in wind_resource_data + has_upperbound = f"{wind_resource_spec}_{height_upper}m" in wind_resource_data + if not has_lowerbound or not has_upperbound: + msg = ( + f"Wind resource data for {wind_resource_spec} is missing either " + f"{height_lower}m or {height_upper}m" + ) + raise ValueError(msg) + + combined_data = np.stack( + [ + wind_resource_data[f"{wind_resource_spec}_{height_lower}m"], + wind_resource_data[f"{wind_resource_spec}_{height_upper}m"], + ] + ) + averaged_data = combined_data.mean(axis=0) + + return averaged_data diff --git a/h2integrate/converters/wind/tools/test/test_resource_tools.py b/h2integrate/converters/wind/tools/test/test_resource_tools.py new file mode 100644 index 000000000..50a7c2dd0 --- /dev/null +++ b/h2integrate/converters/wind/tools/test/test_resource_tools.py @@ -0,0 +1,157 @@ +import numpy as np +import pytest +import CoolProp +import openmdao.api as om +from pytest import fixture + +from h2integrate.converters.wind.tools.resource_tools import ( + calculate_air_density, + average_wind_data_for_hubheight, + weighted_average_wind_data_for_hubheight, +) +from h2integrate.resource.wind.nrel_developer_wtk_api import WTKNRELDeveloperAPIWindResource + + +@fixture +def wind_resource_data(): + plant_config = { + "site": { + "latitude": 34.22, + "longitude": -102.75, + "resources": { + "wind_resource": { + "resource_model": "wind_toolkit_v2_api", + "resource_parameters": { + "latitude": 35.2018863, + "longitude": -101.945027, + "resource_year": 2012, # 2013, + }, + } + }, + }, + "plant": { + "plant_life": 30, + "simulation": { + "dt": 3600, + "n_timesteps": 8760, + "start_time": "01/01/1900 00:30:00", + "timezone": 0, + }, + }, + } + + prob = om.Problem() + comp = WTKNRELDeveloperAPIWindResource( + plant_config=plant_config, + resource_config=plant_config["site"]["resources"]["wind_resource"]["resource_parameters"], + driver_config={}, + ) + prob.model.add_subsystem("resource", comp) + prob.setup() + prob.run_model() + wtk_data = prob.get_val("resource.wind_resource_data") + + return wtk_data + + +def test_air_density_calcs(subtests): + z = 0 + T = 288.15 - 0.0065 * z + P = 101325 * (1 - 2.25577e-5 * z) ** 5.25588 + rho0 = CoolProp.CoolProp.PropsSI("D", "T", T, "P", P, "Air") + rho0_calc = calculate_air_density(z) + with subtests.test("air density at sea level"): + assert pytest.approx(rho0_calc, abs=1e-3) == rho0 + + z = 500 + T = 288.15 - 0.0065 * z + P = 101325 * (1 - 2.25577e-5 * z) ** 5.25588 + rho500m = CoolProp.CoolProp.PropsSI("D", "T", T, "P", P, "Air") + rho500m_calc = calculate_air_density(z) + with subtests.test("air density at 1000m"): + assert pytest.approx(rho500m_calc, abs=1e-3) == rho500m + + +def test_resource_averaging(wind_resource_data, subtests): + avg_windspeed = average_wind_data_for_hubheight(wind_resource_data, [100, 140], "wind_speed") + i_lb_less_than_ub = np.argwhere( + wind_resource_data["wind_speed_100m"] < wind_resource_data["wind_speed_140m"] + ).flatten() + with subtests.test("100m wind speed < average wind speed"): + np.testing.assert_array_less( + wind_resource_data["wind_speed_100m"][i_lb_less_than_ub], + avg_windspeed[i_lb_less_than_ub], + ) + + with subtests.test("140m wind speed > average wind speed"): + np.testing.assert_array_less( + avg_windspeed[i_lb_less_than_ub], + wind_resource_data["wind_speed_140m"][i_lb_less_than_ub], + ) + + with subtests.test("Wind speed at t=0 is average"): + mean_ws = np.mean( + [wind_resource_data["wind_speed_100m"][0], wind_resource_data["wind_speed_140m"][0]] + ) + assert pytest.approx(avg_windspeed[0], rel=1e-6) == mean_ws + + with subtests.test("Avg Wind speed at t=0"): + assert pytest.approx(avg_windspeed[0], rel=1e-6) == 15.78 + + +def test_resource_equal_weighted_averaging(wind_resource_data, subtests): + hub_height = 120 + avg_windspeed = weighted_average_wind_data_for_hubheight( + wind_resource_data, [100, 140], hub_height, "wind_speed" + ) + i_lb_less_than_ub = np.argwhere( + wind_resource_data["wind_speed_100m"] < wind_resource_data["wind_speed_140m"] + ).flatten() + with subtests.test("100m wind speed < average wind speed"): + np.testing.assert_array_less( + wind_resource_data["wind_speed_100m"][i_lb_less_than_ub], + avg_windspeed[i_lb_less_than_ub], + ) + + with subtests.test("140m wind speed > average wind speed"): + np.testing.assert_array_less( + avg_windspeed[i_lb_less_than_ub], + wind_resource_data["wind_speed_140m"][i_lb_less_than_ub], + ) + + with subtests.test("Wind speed at t=0 is average"): + mean_ws = np.mean( + [wind_resource_data["wind_speed_100m"][0], wind_resource_data["wind_speed_140m"][0]] + ) + assert pytest.approx(avg_windspeed[0], rel=1e-6) == mean_ws + + with subtests.test("Avg Wind speed at t=0"): + assert pytest.approx(avg_windspeed[0], rel=1e-6) == 15.78 + + +def test_resource_unequal_weighted_averaging(wind_resource_data, subtests): + hub_height = 135 + weighted_avg_windspeed = weighted_average_wind_data_for_hubheight( + wind_resource_data, [100, 140], hub_height, "wind_speed" + ) + + lb_to_avg_diff = np.abs(weighted_avg_windspeed - wind_resource_data["wind_speed_100m"]) + + avg_to_ub_diff = np.abs(weighted_avg_windspeed - wind_resource_data["wind_speed_140m"]) + + i_nonzero_diff = np.argwhere( + (wind_resource_data["wind_speed_100m"] - wind_resource_data["wind_speed_140m"]) != 0 + ).flatten() + + # the weighted_avg_windspeed should be closer to the 140m height than the 100m + with subtests.test("Weighted avg wind speed is closer to 140m wind speed than 100m wind speed"): + np.testing.assert_array_less(avg_to_ub_diff[i_nonzero_diff], lb_to_avg_diff[i_nonzero_diff]) + + with subtests.test("Weighted avg wind speed at t=0 is greater than equal average"): + mean_ws = np.mean( + [wind_resource_data["wind_speed_100m"][0], wind_resource_data["wind_speed_140m"][0]] + ) + assert weighted_avg_windspeed[0] > mean_ws + + with subtests.test("Avg Wind speed at t=0"): + assert pytest.approx(weighted_avg_windspeed[0], rel=1e-6) == 16.56 diff --git a/h2integrate/core/model_baseclasses.py b/h2integrate/core/model_baseclasses.py index e9ea6c4d9..296cc214a 100644 --- a/h2integrate/core/model_baseclasses.py +++ b/h2integrate/core/model_baseclasses.py @@ -164,7 +164,7 @@ def __attrs_post_init__(self): self.cache_dir = Path(self.cache_dir) # Create a cache directory if it doesn't exist - if not self.cache_dir.exists(): + if self.enable_caching and not self.cache_dir.exists(): self.cache_dir.mkdir(parents=True, exist_ok=True) diff --git a/h2integrate/core/supported_models.py b/h2integrate/core/supported_models.py index c9b41974a..466c590c3 100644 --- a/h2integrate/core/supported_models.py +++ b/h2integrate/core/supported_models.py @@ -6,6 +6,7 @@ from h2integrate.finances.profast_lco import ProFastLCO from h2integrate.finances.profast_npv import ProFastNPV from h2integrate.converters.steel.steel import SteelPerformanceModel, SteelCostAndFinancialModel +from h2integrate.converters.wind.floris import FlorisWindPlantPerformanceModel from h2integrate.converters.iron.iron_mine import ( IronMineCostComponent, IronMinePerformanceComponent, @@ -170,6 +171,7 @@ # Converters "atb_wind_cost": ATBWindPlantCostModel, "pysam_wind_plant_performance": PYSAMWindPlantPerformanceModel, + "floris_wind_plant_performance": FlorisWindPlantPerformanceModel, "pysam_solar_plant_performance": PYSAMSolarPlantPerformanceModel, "atb_utility_pv_cost": ATBUtilityPVCostModel, "atb_comm_res_pv_cost": ATBResComPVCostModel, diff --git a/h2integrate/preprocess/test/test_wind_turbine_file_tools.py b/h2integrate/preprocess/test/test_wind_turbine_file_tools.py index 772010619..21735265a 100644 --- a/h2integrate/preprocess/test/test_wind_turbine_file_tools.py +++ b/h2integrate/preprocess/test/test_wind_turbine_file_tools.py @@ -4,6 +4,7 @@ from h2integrate import EXAMPLE_DIR from h2integrate.core.utilities import load_yaml +from h2integrate.converters.wind.floris import FlorisWindPlantPerformanceModel from h2integrate.core.inputs.validation import load_tech_yaml, load_plant_yaml from h2integrate.converters.wind.wind_pysam import PYSAMWindPlantPerformanceModel from h2integrate.preprocess.wind_turbine_file_tools import ( @@ -97,6 +98,65 @@ def test_floris_turbine_export(subtests): assert output_fpath.exists() assert output_fpath.is_file() - # TODO: add test with floris model - # with subtests.test("File runs with FLORIS model"): - # pass + floris_options = load_yaml(output_fpath) + + plant_config_path = EXAMPLE_DIR / "05_wind_h2_opt" / "plant_config.yaml" + tech_config_path = EXAMPLE_DIR / "26_floris" / "tech_config.yaml" + + plant_config = load_plant_yaml(plant_config_path) + tech_config = load_tech_yaml(tech_config_path) + + updated_parameters = { + "hub_height": -1, + "floris_turbine_config": floris_options, + } + + tech_config["technologies"]["wind"]["model_inputs"]["performance_parameters"].update( + updated_parameters + ) + + prob = om.Problem() + wind_resource = WTKNRELDeveloperAPIWindResource( + plant_config=plant_config, + resource_config=plant_config["site"]["resources"]["wind_resource"]["resource_parameters"], + driver_config={}, + ) + + wind_plant = FlorisWindPlantPerformanceModel( + plant_config=plant_config, + tech_config=tech_config["technologies"]["wind"], + driver_config={}, + ) + + prob.model.add_subsystem("wind_resource", wind_resource, promotes=["*"]) + prob.model.add_subsystem("wind_plant", wind_plant, promotes=["*"]) + prob.setup() + prob.run_model() + + with subtests.test("File runs with Floris, check total capacity"): + assert ( + pytest.approx(prob.get_val("wind_plant.total_capacity", units="MW"), rel=1e-6) == 600.0 + ) + + with subtests.test("File runs with Floris, check turbine size"): + assert ( + pytest.approx(prob.get_val("wind_plant.num_turbines", units="unitless"), rel=1e-6) + == 100 + ) + + with subtests.test("File runs with Floris, check hub-height"): + assert pytest.approx(prob.get_val("wind_plant.hub_height", units="m"), rel=1e-6) == 140.0 + + with subtests.test("File runs with Floris, check capacity factor"): + assert ( + pytest.approx(prob.get_val("wind_plant.capacity_factor", units="percent")[0], rel=1e-6) + == 53.556784 + ) + + with subtests.test("File runs with Floris, check AEP"): + assert ( + pytest.approx( + prob.get_val("wind_plant.total_electricity_produced", units="MW*h/yr")[0], rel=1e-6 + ) + == 2814944.574 + ) diff --git a/h2integrate/resource/wind/openmeteo_wind.py b/h2integrate/resource/wind/openmeteo_wind.py index f8b5edac0..61574b37e 100644 --- a/h2integrate/resource/wind/openmeteo_wind.py +++ b/h2integrate/resource/wind/openmeteo_wind.py @@ -80,6 +80,7 @@ def setup(self): "surface_pressure": "hPa", "precipitation": "mm/h", "relative_humidity_2m": "unitless", + "is_day": "percent", } # get the data dictionary data = self.get_data(self.config.latitude, self.config.longitude) @@ -319,6 +320,10 @@ def format_timeseries_data(self, data): if old_c not in self.hourly_wind_data_to_units: continue + if "is_day" in c: + data_rename_mapper.update({c: "is_day"}) + data_units.update({"is_day": "percent"}) + if "surface" in c: new_c += "_0m" new_c = new_c.replace("surface", "").replace("__", "").strip("_") diff --git a/h2integrate/resource/wind/wind_resource_base.py b/h2integrate/resource/wind/wind_resource_base.py index cd6a7f766..f45876a97 100644 --- a/h2integrate/resource/wind/wind_resource_base.py +++ b/h2integrate/resource/wind/wind_resource_base.py @@ -14,6 +14,7 @@ def setup(self): "pressure": "atm", "precipitation_rate": "mm/h", "relative_humidity": "percent", + "is_day": "unitless", } def compare_units_and_correct(self, data, data_units): diff --git a/library/floris_turbine_Vestas_660kW.yaml b/library/floris_turbine_Vestas_660kW.yaml new file mode 100644 index 000000000..e7fa46a8f --- /dev/null +++ b/library/floris_turbine_Vestas_660kW.yaml @@ -0,0 +1,213 @@ +turbine_type: VestasV47_660kW_47 +hub_height: 65.0 +TSR: 8.0 +rotor_diameter: 47.0 +power_thrust_table: + ref_air_density: 1.225 + ref_tilt: 5.0 + cosine_loss_exponent_yaw: 1.88 + cosine_loss_exponent_tilt: 1.88 + wind_speed: + - 0.0 + - 1.0 + - 2.0 + - 3.0 + - 4.0 + - 4.17 + - 4.58 + - 5.02 + - 5.5 + - 6.0 + - 6.51 + - 6.99 + - 7.49 + - 7.99 + - 8.47 + - 9.02 + - 9.51 + - 10.01 + - 10.47 + - 10.99 + - 11.5 + - 11.99 + - 12.48 + - 13.01 + - 13.49 + - 13.99 + - 14.5 + - 14.91 + - 15.39 + - 15.96 + - 16.42 + - 16.91 + - 17.45 + - 17.91 + - 17.91 + - 18.91 + - 19.91 + - 20.91 + - 21.91 + - 22.91 + - 23.91 + - 24.91 + - 25.91 + - 26.91 + - 27.91 + - 28.91 + - 29.91 + - 30.91 + - 31.91 + - 32.91 + - 33.91 + - 34.91 + - 35.91 + - 36.91 + - 37.91 + - 38.91 + - 39.91 + - 40.91 + - 41.91 + - 42.91 + - 43.91 + - 44.91 + - 45.91 + - 46.91 + - 47.91 + - 48.91 + - 49.91 + power: + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 11.75 + - 26.46 + - 41.79 + - 67.86 + - 95.69 + - 131.45 + - 164.31 + - 209.85 + - 253.8 + - 301.19 + - 362.48 + - 406.7 + - 468.74 + - 508.55 + - 551.26 + - 587.82 + - 612.36 + - 634.38 + - 646.79 + - 651.4 + - 658.4 + - 659.46 + - 660.0 + - 660.0 + - 660.0 + - 660.0 + - 660.0 + - 660.0 + - 660.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + thrust_coefficient: + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.1596 + - 0.2755 + - 0.3399 + - 0.4364 + - 0.4857 + - 0.529 + - 0.5382 + - 0.5653 + - 0.5637 + - 0.5604 + - 0.5572 + - 0.5274 + - 0.5183 + - 0.4857 + - 0.4488 + - 0.412 + - 0.3728 + - 0.3386 + - 0.3018 + - 0.2696 + - 0.2416 + - 0.2164 + - 0.1984 + - 0.1794 + - 0.1596 + - 0.1466 + - 0.1336 + - 0.1208 + - 0.1112 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 + - 0.0 diff --git a/library/floris_v4_default_template.yaml b/library/floris_v4_default_template.yaml new file mode 100644 index 000000000..4797b6083 --- /dev/null +++ b/library/floris_v4_default_template.yaml @@ -0,0 +1,91 @@ + +name: Gauss +description: Onshore template +floris_version: v4.0.0 +logging: + console: + enable: false + level: WARNING + file: + enable: false + level: WARNING +solver: + type: turbine_grid + turbine_grid_points: 1 +flow_field: + air_density: 1.225 + reference_wind_height: -1 + wind_directions: [] + wind_shear: 0.33 + wind_speeds: [] + wind_veer: 0.0 + +wake: + model_strings: + combination_model: sosfs + deflection_model: gauss + turbulence_model: crespo_hernandez + velocity_model: gauss + enable_secondary_steering: false + enable_yaw_added_recovery: false + enable_transverse_velocities: false + wake_deflection_parameters: + gauss: + ad: 0.0 + alpha: 0.58 + bd: 0.0 + beta: 0.077 + dm: 1.0 + ka: 0.38 + kb: 0.004 + jimenez: + ad: 0.0 + bd: 0.0 + kd: 0.05 + wake_velocity_parameters: + cc: + a_s: 0.179367259 + b_s: 0.0118889215 + c_s1: 0.0563691592 + c_s2: 0.13290157 + a_f: 3.11 + b_f: -0.68 + c_f: 2.41 + alpha_mod: 1.0 + gauss: + alpha: 0.58 + beta: 0.077 + ka: 0.38 + kb: 0.004 + jensen: + we: 0.05 + wake_turbulence_parameters: + crespo_hernandez: + initial: 0.1 + constant: 0.5 + ai: 0.8 + downstream: -0.32 + enable_active_wake_mixing: false + + wake_velocity_parameters: + cc: + a_f: 3.11 + a_s: 0.179367259 + alpha_mod: 1.0 + b_f: -0.68 + b_s: 0.0118889215 + c_f: 2.41 + c_s1: 0.0563691592 + c_s2: 0.13290157 + gauss: + alpha: 0.58 + beta: 0.077 + ka: 0.38 + kb: 0.004 + jensen: + we: 0.05 +farm: + layout_x: [] + layout_y: [] + # turbine_type: + # - operation_model: cosine-loss