guanlan-src/scripts/sim_hub/wake.py

"""尾流示意: 同一条山脊, 换大叶轮后排布逻辑整个变了。
单一风向 270°、机位完全对齐的示意布局(非任何真实风场)。
引擎 = PyWake (DTU, MIT) · Bastankhah-Porte-Agel 2014 高斯尾流模型。"""
import os, numpy as np, matplotlib
matplotlib.use("Agg")
import matplotlib.pyplot as plt
from matplotlib import font_manager as fm
from py_wake.site import UniformSite
from py_wake.wind_turbines import WindTurbine
from py_wake.wind_turbines.power_ct_functions import PowerCtTabular
from py_wake.literature.gaussian_models import Bastankhah_PorteAgel_2014

for cand in ["PingFang SC", "Hiragino Sans GB", "Songti SC", "STHeiti"]:
    try:
        fm.findfont(fm.FontProperties(family=cand), fallback_to_default=False)
        matplotlib.rcParams["font.sans-serif"] = [cand, "DejaVu Sans"]; break
    except Exception: continue
matplotlib.rcParams["axes.unicode_minus"] = False

BG="#04080E"; INK="#EAF1F7"; INK2="#A9BECD"; MUT="#7C93A6"; AMB="#F0C674"; GRN="#7FD9B8"
RIDGE = 6300.   # 一条 6.3 km 的山脊

def mk_wt(name, D, H, rated_kw):
    u = np.arange(0, 26.)
    p = np.clip(((u-3)/(11-3))**3, 0, 1)*rated_kw*1e3; p[u<3]=0; p[u>25]=0
    ct = np.where((u>=3)&(u<=25), np.clip(0.8*np.minimum(1,(11/np.maximum(u,3))**2), .05, .8), 0)
    return WindTurbine(name=name, diameter=D, hub_height=H,
                       powerCtFunction=PowerCtTabular(u, p, 'w', ct))

site = UniformSite(p_wd=[1], ti=.10, ws=9.)
wt15 = mk_wt("1.5MW", 82., 70., 1500.)
wt50 = mk_wt("5MW", 195., 110., 5550.)

cases = [
    ("①", "1.5 MW · 叶轮 82 m",  wt15, 82.,  500.,  GRN, "间距 500 m ＝ 6.1 D　尾流有恢复余地"),
    ("②", "5 MW · 叶轮 195 m",   wt50, 195., 500.,  AMB, "间距 500 m ＝ 2.6 D　同样的米数，完全不够"),
    ("③", "5 MW · 叶轮 195 m",   wt50, 195., 1190., GRN, "间距 1190 m ＝ 6.1 D　恢复到与 ① 同等条件"),
]

fig = plt.figure(figsize=(14.6, 7.6), facecolor=BG)
gs = fig.add_gridspec(3, 1, left=.052, right=.895, top=.845, bottom=.085, hspace=.62)

for i, (num, name, wt, D, sp, col, note) in enumerate(cases):
    ax = fig.add_subplot(gs[i]); ax.set_facecolor(BG)
    n = int(RIDGE // sp) + 1
    x = np.arange(n)*sp; y = np.zeros(n)
    wfm = Bastankhah_PorteAgel_2014(site, wt, k=0.0324555)
    sim = wfm(x, y, wd=270, ws=9.)
    fmap = sim.flow_map(wd=270, ws=9.)
    Z = fmap.WS_eff.squeeze().values; X = fmap.x.values; Y = fmap.y.values
    if Z.shape[0] == len(X): Z = Z.T
    c = ax.contourf(X, Y, Z, levels=np.linspace(5.0, 9.2, 64), cmap="magma", extend="both")
    for xi in x:
        ax.plot([xi, xi], [-D/2, D/2], color=INK, lw=2.0, solid_capstyle="round", zorder=5)
        ax.plot(xi, 0, marker="o", ms=3.2, color=INK, zorder=6)

    # 单行标题, 顺序拼接, 不会重叠
    ax.text(0, 1.44, num, transform=ax.transAxes, color=col, fontsize=15, weight="bold", va="top")
    ax.text(.022, 1.44, name, transform=ax.transAxes, color=INK, fontsize=13.5, weight="bold", va="top")
    ax.text(.022, 1.44, "　" + " "*len(name)*2, transform=ax.transAxes, alpha=0, fontsize=13.5, va="top")
    ax.text(.215, 1.44, note, transform=ax.transAxes, color=col, fontsize=12.5, va="top")
    ax.text(.86, 1.44, f"这条山脊放得下 {n} 台", transform=ax.transAxes,
            color=INK if i != 2 else GRN, fontsize=13, weight="bold", va="top", ha="left")

    ax.set_xlim(-320, RIDGE + 320); ax.set_ylim(-300, 300); ax.set_aspect("equal")
    for s in ax.spines.values(): s.set_color("#1B2A36")
    ax.tick_params(colors=MUT, labelsize=8.5)
    if i < 2: ax.set_xticklabels([])
    ax.set_ylabel("横向 m", color=MUT, fontsize=9)
ax.set_xlabel("沿山脊顺风向 (m) — 来流 9.0 m/s · 湍流强度 10%", color=MUT, fontsize=9.5)

cax = fig.add_axes([.912, .30, .011, .40])
cb = fig.colorbar(c, cax=cax)
cb.set_label("尾流后的有效风速 (m/s)", color=INK2, fontsize=9.5)
cb.ax.yaxis.set_tick_params(color=MUT, labelsize=8.5)
plt.setp(plt.getp(cb.ax.axes, 'yticklabels'), color=MUT)
cb.outline.set_edgecolor("#1B2A36")

fig.text(.052, .962, "换大叶轮不是换个大的就行 —— 整条山脊的排布逻辑都变了",
         color=INK, fontsize=19, weight="bold")
fig.text(.052, .916, "同一条 6.3 km 山脊。要维持同等尾流条件，5 MW 机组的间距得拉到 1190 m，机位数从 13 台掉到 6 台。",
         color=INK2, fontsize=12)
fig.text(.052, .888, "单一风向、机位完全对齐的示意布局（非任何真实风场，推力曲线为通用示意值）；实际风场有风向分布，尾流损失远小于本图，此处只为说明原理。",
         color=MUT, fontsize=9.3)
fig.text(.052, .022, "引擎：PyWake（丹麦技术大学，MIT 许可）· Bastankhah–Porté-Agel 2014 高斯尾流模型",
         color="#5E7385", fontsize=8.8)

out = os.path.join(os.path.dirname(__file__), "wake_compare.png")
fig.savefig(out, dpi=150, facecolor=BG)
print("SAVED", out)