Python Wind Farm Analysis — Deep Dive

Technical foundation

Wind farm analysis in Python spans fluid dynamics, statistics, optimization, and time-series forecasting. Production-grade analyses must handle mesoscale weather models, complex terrain effects, multi-wake superposition, and uncertainty quantification. This deep dive covers the computational techniques and code patterns used by wind energy professionals.

Wind data processing and Weibull fitting

Raw met mast or LiDAR data arrives as 10-minute averages of speed, direction, temperature, and pressure at multiple heights. The first task is quality control and statistical characterization.

import numpy as np
import pandas as pd
from scipy.stats import weibull_min
from scipy.optimize import minimize_scalar

def fit_weibull(speeds: np.ndarray) -> tuple[float, float]:
    """Fit Weibull distribution using maximum likelihood."""
    speeds = speeds[speeds > 0]  # Exclude calm periods
    shape, loc, scale = weibull_min.fit(speeds, floc=0)
    return shape, scale

def wind_power_density(speeds: np.ndarray, air_density: float = 1.225) -> float:
    """Calculate mean wind power density in W/m²."""
    return 0.5 * air_density * np.mean(speeds ** 3)

# Load and process met mast data
mast = pd.read_csv("met_mast_80m.csv", parse_dates=["timestamp"], index_col="timestamp")
mast = mast[mast["speed_80m"] > 0].dropna(subset=["speed_80m", "direction_80m"])

k, A = fit_weibull(mast["speed_80m"].values)
wpd = wind_power_density(mast["speed_80m"].values)
print(f"Weibull k={k:.2f}, A={A:.2f} m/s, WPD={wpd:.0f} W/m²")

Wind power density above 400 W/m² at hub height is generally considered excellent for commercial development.

Vertical wind shear extrapolation

Measurements at one height need extrapolation to hub height. The power law is standard:

def extrapolate_wind_speed(
    speed_measured: np.ndarray,
    height_measured: float,
    height_target: float,
    shear_exponent: float = 0.143,  # 1/7 power law default
) -> np.ndarray:
    """Extrapolate wind speed to target height using power law."""
    return speed_measured * (height_target / height_measured) ** shear_exponent

def calculate_shear_exponent(
    speed_lower: np.ndarray, height_lower: float,
    speed_upper: np.ndarray, height_upper: float,
) -> float:
    """Calculate site-specific shear from two measurement heights."""
    mask = (speed_lower > 1) & (speed_upper > 1)
    ratio = np.log(speed_upper[mask] / speed_lower[mask])
    height_ratio = np.log(height_upper / height_lower)
    return np.median(ratio / height_ratio)

Typical onshore shear exponents range from 0.10 (flat, open terrain) to 0.30 (forested or complex terrain). Using the default 1/7 rule when the site has high shear can underestimate hub-height wind speed by 10–15%.

Measure-Correlate-Predict (MCP)

Short on-site campaigns (1–2 years) must be extended to represent long-term conditions (10–20 years) using nearby reference data:

from sklearn.linear_model import LinearRegression

def mcp_linear(site_data: pd.Series, reference_data: pd.Series,
               reference_longterm: pd.Series) -> pd.Series:
    """Simple linear MCP to predict long-term site wind from reference."""
    # Concurrent period
    concurrent = pd.concat([site_data, reference_data], axis=1).dropna()
    concurrent.columns = ["site", "ref"]

    model = LinearRegression()
    model.fit(concurrent[["ref"]], concurrent["site"])

    # Predict long-term
    ref_lt = reference_longterm.to_frame(name="ref")
    predicted = model.predict(ref_lt)
    return pd.Series(predicted.flatten(), index=reference_longterm.index)

More sophisticated MCP methods (variance ratio, matrix method, Weibull-based) better preserve the site’s wind speed distribution. The brightwind library implements several MCP approaches.

Wake modeling with FLORIS

NREL’s FLORIS is the most actively developed Python wake modeling framework:

import floris

# Initialize FLORIS with Gauss-Curl-Hybrid wake model
fi = floris.FlorisInterface("gch.yaml")

# Define a 3x3 grid layout
x_positions = [0, 0, 0, 630, 630, 630, 1260, 1260, 1260]
y_positions = [0, 630, 1260, 0, 630, 1260, 0, 630, 1260]

fi.reinitialize(
    layout_x=x_positions,
    layout_y=y_positions,
    wind_speeds=[8.0, 10.0, 12.0],
    wind_directions=[270.0, 280.0, 290.0],
)

fi.calculate_wake()
turbine_powers = fi.get_turbine_powers()  # Watts per turbine per condition
farm_power = fi.get_farm_power()  # Total farm power per condition

Multi-wake superposition

When a downstream turbine sits in multiple overlapping wakes, superposition models determine the combined effect. Common approaches:

  • Linear (Lissaman) — Velocity deficits add linearly. Simple but overestimates losses.
  • Sum of squares (Katic)deficit_combined = sqrt(sum(deficit_i²)). The most widely used method.
  • Max deficit — Takes the largest single wake. Underestimates combined effects.

FLORIS supports all three and allows custom superposition functions.

Layout optimization

FLORIS includes optimization tools for layout and yaw angle:

from floris.tools.optimization.yaw_optimization.yaw_optimizer_sr import YawOptimizationSR

# Yaw optimization — angle turbines slightly off-wind to deflect wakes
yaw_opt = YawOptimizationSR(fi)
df_opt = yaw_opt.optimize()
optimal_yaw_angles = df_opt["yaw_angles_opt"].values

# Layout optimization using scipy
from floris.tools.optimization.layout_optimization.layout_optimization_scipy import (
    LayoutOptimizationScipy,
)

# Define site boundary as polygon
boundary = [(0, 0), (3000, 0), (3000, 3000), (0, 3000)]

layout_opt = LayoutOptimizationScipy(
    fi,
    boundaries=boundary,
    min_dist=400,  # 4D minimum spacing for 100m rotor
)
sol = layout_opt.optimize()

Yaw-based wake steering

Modern turbines can intentionally yaw (turn) slightly off the incoming wind direction. This reduces the yawed turbine’s output by 1–3% but deflects its wake away from downstream turbines, increasing their output by 5–15%. The net farm-level gain is typically 1–4% AEP. This is an active area of research where FLORIS is heavily used.

Annual energy production calculation

AEP combines wind frequency distribution with wake-affected power output:

def calculate_aep(
    wind_speeds: np.ndarray,
    wind_directions: np.ndarray,
    frequency_table: np.ndarray,
    floris_interface,
    hours_per_year: float = 8760,
) -> float:
    """Calculate AEP from wind rose frequency table and FLORIS."""
    floris_interface.reinitialize(
        wind_speeds=wind_speeds,
        wind_directions=wind_directions,
    )
    floris_interface.calculate_wake()
    farm_powers = floris_interface.get_farm_power()  # W per (direction, speed) bin

    # Weight by frequency and convert to MWh
    aep_mwh = np.sum(farm_powers * frequency_table) * hours_per_year / 1e6
    return aep_mwh

Uncertainty quantification

AEP estimates carry significant uncertainty. Industry practice quantifies P50 (median), P75, and P90 estimates:

def aep_uncertainty(
    aep_central: float,
    wind_resource_uncertainty: float = 0.06,  # 6% typical
    wake_model_uncertainty: float = 0.05,
    power_curve_uncertainty: float = 0.03,
    availability_uncertainty: float = 0.02,
) -> dict:
    """Calculate P50, P75, P90 AEP estimates."""
    # Combine independent uncertainties in quadrature
    total_sigma = np.sqrt(
        wind_resource_uncertainty**2
        + wake_model_uncertainty**2
        + power_curve_uncertainty**2
        + availability_uncertainty**2
    )
    from scipy.stats import norm
    return {
        "P50": aep_central,
        "P75": aep_central * (1 - norm.ppf(0.75) * total_sigma),
        "P90": aep_central * (1 - norm.ppf(0.90) * total_sigma),
        "total_uncertainty": total_sigma,
    }

Banks typically finance projects based on P75 or P90 estimates, not P50. A 10% total uncertainty means the P90 estimate is roughly 13% below P50 — a significant gap that affects project economics.

Production forecasting for operations

Once a farm is operating, day-ahead production forecasts support grid dispatch and power trading:

import lightgbm as lgb

def train_production_forecast(historical: pd.DataFrame):
    """Train short-term production forecast model."""
    features = [
        "wind_speed_forecast", "wind_direction_forecast",
        "temperature_forecast", "pressure_forecast",
        "hour", "month", "production_lag_1h", "production_lag_24h",
    ]
    target = "farm_production_mw"

    model = lgb.LGBMRegressor(n_estimators=500, learning_rate=0.05)
    model.fit(historical[features], historical[target])
    return model

Combining NWP (Numerical Weather Prediction) forecasts with autoregressive features and SCADA data typically achieves 5–10% MAPE for 1–24 hour horizons.

Tradeoffs

DecisionOption AOption B
Wake modelJensen (fast, conservative)Gauss-Curl-Hybrid (accurate, slower)
Layout approachRegular grid (predictable)Optimized layout (higher AEP, complex permitting)
Measurement campaign1 year (cheaper)2+ years (captures inter-annual variability)
Turbine spacing4D (high density, more wake loss)8D+ (low density, less land utilization)
Yaw controlGreedy (simple)Coordinated wake steering (1–4% gain, requires advanced controls)

One thing to remember: Wind farm analysis is an optimization problem where uncertainty quantification is as important as the point estimate — banks, insurers, and grid operators all make decisions based on P75/P90 numbers, not best-case scenarios.

pythonwind-energydata-sciencesustainability

See Also