改进奇诺多面体的需求侧资源可行域聚合方法
2026/9/19 15:28:37 网站建设 项目流程

简介:资源为基于改进奇诺多面体的需求侧资源可行域聚合研究论文复现包,适合电力系统研究人员、需求侧管理工程师以及对优化建模和调度算法感兴趣的读者。核心解决柔性负荷、储能、电动汽车等分散资源容量小、特性差异大、难以统一聚合的难题,内容涵盖奇诺多面体基本概念、生成器改进设计(如储能损耗、充电时间约束、功率限制)、闵可夫斯基和聚合方法、PCA降维处理维数灾难以及多方向优化提高近似精度等。资源共1个文件,为PDF文档,整体约793KB,内含论文原文及配套Python代码复现与详细解释,代码中实现了ImprovedZonotope类,包含采样、绘图、闵可夫斯基和、改进生成器等功能模块,便于直接运行和二次开发。目前已有176人学习浏览,是电力系统需求侧资源聚合建模与工程应用方面较为实用的参考资料。

1. 用奇诺多面体做需求侧可行域聚合,开局先接受“不精确”

配网调度员每天要面对上百个储能站、空调集群和充电桩的联合调节问题。每个设备的可调度集合是一堆线性不等式围出的高维多面体,把这些集合做闵可夫斯基和,多面体的顶点数以组合速度膨胀,直接交给求解器计算会在预处理阶段就卡住。奇诺多面体(zonotope)用“中心 + 生成器矩阵”编码凸体,让求和退化成矩阵拼接,计算量从顶点枚举降到矩阵乘法级别;论文里针对它做的“改进”,通常落在生成器缩减、有界约束、主方向挑选上,把聚合结果压到几十个生成器仍保持误差可控。这篇文章围绕“基于改进奇诺多面体的需求侧资源可行域聚合研究”展开,先立数学模型,再给一套可在 Python 中运行的复现代码,覆盖资源建模、闵可夫斯基和聚合、生成器缩减、二维投影与误差验证。适合要复现同类论文却拿不到开源代码的工程师,也适合想在自有负荷数据上跑通可行域聚合流程的调度侧开发者。

2. 需求侧可行域的数学模型与奇诺多面体聚合封闭性

2.1 从设备约束到可调度集合:可行域聚合在算什么东西

需求侧资源里最常见的一类模型是储能,它的约束由功率限和能量限共同组成。以一台额定功率 P_max、容量 E_max、初始电量为 E0 的电池为例,调度间隔取 dt,T 个时段内的有功功率序列 p ∈ R^T,约束写成:

import numpy as np # 构造单台储能的约束矩阵 A x <= b # x = [p_0, p_1, ..., p_{T-1}],单位 MW T = 96 # 调度时段数 dt = 0.25 # 小时,15 分钟间隔 P_max = 1.0 # MW E_max = 4.0 # MWh E0 = 1.5 # MWh A_p = np.vstack([np.eye(T), -np.eye(T)]) b_p = np.r_[np.full(T, P_max), np.full(T, P_max)] L = np.tril(np.ones((T, T))) * dt # 电量累积矩阵:E0 + L @ p A_e = np.vstack([L, -L]) b_e = np.r_[np.full(T, E_max - E0), np.full(T, E0)] A = np.vstack([A_p, A_e]) b = np.hstack([b_p, b_e])

A 的行数在 T=96 时达到 288 行,每个资源的可行域是 R^96 里的一个高维多面体。多个资源聚合时,真正要算的是闵可夫斯基和:Θ_agg = Θ_1 ⊕ Θ_2 ⊕ … ⊕ Θ_N = { Σ_i x^(i) : x^(i) ∈ Θ_i }。这个集合的意义是,调度中心只需要向聚合商下发一条总功率曲线,聚合商能把它分解到每台设备上执行。问题在于多面体的闵可夫斯基和需要顶点枚举,单台的顶点数已随时段数增长,N 台聚合后几乎不可解。

2.2 奇诺多面体的定义与三条封闭运算

奇诺多面体是一种特殊的凸多面体,定义为向量 c 与一组一维线段做闵可夫斯基和的结果:

Z = { c + G z : ‖z‖∞ ≤ 1 }

其中 c ∈ R^d 是中心,G ∈ R^(d×m) 的每一列都是一个生成器,z 被限制在无限范数球内。直观理解就是:从中心出发,沿着 m 个方向各伸出一段长度,这些线段的所有组合叠加出的凸体。它的关键优势在于封闭性,下面三个性质是聚合算法全部依赖的数学基础:

运算公式复杂度
线性映射R·Z = { Rc + RG z }O(rdm)
闵可夫斯基和Z_1 ⊕ Z_2 = { (c_1+c_2) + [G_1, G_2] z }O(d) 拼接
笛卡尔积Z_1 × Z_2,生成器按块对角放置O(d) 拼接

代码里对应一个很短的类实现:

class Zonotope: """中心 c: (d,1),生成器矩阵 G: (d,m)""" def __init__(self, c, G): self.c = np.asarray(c, dtype=float).reshape(-1, 1) self.G = np.asarray(G, dtype=float) def minkowski_sum(self, other): # 闵可夫斯基和:中心相加,生成器横向拼接 return Zonotope(self.c + other.c, np.hstack([self.G, other.G])) def linear_map(self, R): # 线性变换:左乘一个映射矩阵,常见于投影 return Zonotope(R @ self.c, R @ self.G)

参数上最关键的是生成器矩阵的维度。c 的维度 d 是状态维度,在可行域聚合里通常是时段数;m 是生成器个数,决定了集合的表达精细度。minkowski_sum 不做任何线性规划,只是把两个资源的生成器拼在一起,这正是奇诺多面体能被规模化应用的原因。

2.3 生成器拼接的陷阱:为什么论文都要加“改进”

直接拼接看起来解决了复杂度问题,但代价是生成器数量线性累加。上面每台储能如果近似用 98 个生成器表示,10 台聚合后就是 980 个,100 台就是 9800 个。生成器数量膨胀后,后续再做二次调度优化时,变量 z 的维度变成生成器总数,求解器仍然会慢下来。

更隐蔽的问题出现在表达精度上。单纯拼接只是把所有资源的几何信息“堆”在一起,并没有剔除相互冗余的方向。储能资源的可行域在时段维度上有强关联:相邻时段的功率生成器高度相关,几十个生成器实际上只刻画了两三个独立的方向(充、放、能量偏移)。所以论文里所谓“改进奇诺多面体”,核心动作就是两步:一是有界化,让奇诺多面体表达非对称可行域;二是生成器缩减,在保持集合方向特征的前提下把 m 从上千降到个位数。

3. 改进奇诺多面体:有界参数化与生成器缩减算法

3.1 对称性限制:纯奇诺多面体表达不了一个储能站

按第 2 节公式定义的奇诺多面体关于中心点严格对称。这个性质在数学上很漂亮,但在需求侧资源建模里是硬伤。电池的可行域中,初始电量 E0 通常不在容量区间中间,允许充电的空间和允许放电的空间不对称;空调负荷的舒适度边界、电动汽车的到达离网时间也天然不对称。

改进方向之一是有界奇诺多面体(constrained zonotope),在 z 上额外加一组线性不等式:

Z_c = { c + G z : ‖z‖∞ ≤ 1, L z ≤ d }

L 和 d 是新增的约束矩阵,它们破坏了无限范数球的对称性,让集合可以偏向某个方向。聚合时,中心与生成器的处理方式不变,但 L、d 需要追加拼接:L_agg = blockdiag(L_1, …, L_N),d_agg = [d_1; …; d_N]。这一步在大多数论文的实现里都是直接拼,真正的计算难点不在约束本身,而是约束会让生成器缩减时的包含判定变得复杂。

3.2 外部近似还是内部近似:先把保守方向选对

做生成器缩减之前,要先回答一个问题:聚合后的集合允许比真实可调度集合大,还是必须比它小。这两个选择对应完全不同的工程语义:

模式集合关系调度含义典型适用
inner(内部近似)聚合集 ⊆ 真实集每条指令均可分解到设备,灵活性偏保守现货市场申报、紧急调频
outer(外部近似)聚合集 ⊇ 真实集灵活性不丢,但可能给出不可行指令规划阶段、容量评估

配网侧做实时调度时通常选 inner,因为下发下去的总功率曲线如果分解不了,会直接导致设备越限。论文复现时,这个模式一般作为参数暴露给使用者,代码里用一个字符串变量控制,后续缩减的缩放系数依它取值。

3.3 生成器缩减:SVD 主方向挑选与投影区间重构

生成器缩减的工程化实现有很多种,我在复现中最常用的是“奇异值分解选方向 + 区间盒重构”的启发式方法。它的假设是:聚合后生成器矩阵 G 的奇异值快速衰减,前 k 个奇异向量已经张成了集合的主要方向,剩余方向对形状的贡献很小,可以去掉或者用很小的补偿项吸收。

算法步骤如下:

  1. 对聚合后的 G 做奇异值分解 G = U S V^T。
  2. 计算奇异值平方的累计占比,保留达到阈值 energy_ratio 的前 k 个方向。
  3. 把原生成器投影到主方向张成的子空间:G_proj = U_k^T G,计算每个主方向上的最大绝对投影值。
  4. 构造新的生成器矩阵 G_new = U_k diag(w_1, …, w_k),其中 w_j 取第 j 个主方向的投影半长。
  5. 根据 inner/outer 模式,对 w_j 统一乘一个缩放系数,inner 取 0.9~0.95,outer 取 1.05~1.1。

对应代码:

def reduce_zonotope(Z, k=10, energy_ratio=0.95, mode='inner'): """生成器缩减:SVD 主方向 + 区间投影半长重构""" c, G = Z.c, Z.G d, m = G.shape U, s, _ = np.linalg.svd(G, full_matrices=False) ratio = np.cumsum(s ** 2) / np.sum(s ** 2) k = min(k, m, int(np.searchsorted(ratio, energy_ratio) + 1)) U_k = U[:, :k] # (d, k) 主方向 G_proj = U_k.T @ G # (k, m) 各生成器在主方向上的坐标 w = np.abs(G_proj).max(axis=1) # 每个主方向上的投影半长 w = w * (0.95 if mode == 'inner' else 1.05) return Zonotope(c, U_k @ (w * np.eye(k)))

参数上,k 和 energy_ratio 是相互制约的:k 限制生成器个数上限,energy_ratio 决定实际保留多少个方向。searchsorted 那段是在累计能量曲线上找第一个超过阈值的下标,保证代码在奇异值衰减快时自动少留方向,衰减慢时多留。inner 模式的 0.95 是抵消 SVD 截断误差的工程系数,不保证严格包含关系,严格包含需要用线性规划逐方向校验。

4. 论文复现:基于改进奇诺多面体的可行域聚合完整代码

4.1 实验场景与数据构造

复现实验用 10 台参数不同的储能组成资源集群,调度周期 96 个时段,每个资源的多面体约束由功率上限、容量上限和初始电量决定。对照基准是“真实聚合域”:因为真实聚合的精度无法直接计算,我用支撑函数法求它在多个方向上的边界点,再重建二维投影多边形。

支撑函数的思路是:对任意方向向量 θ,真实聚合域的支撑值为 h(θ) = Σ_i max_{x∈Θ_i} θ^T x,每个最大值用线性规划求;而奇诺多面体的支撑值有解析公式 h(θ) = θ^T c + ‖G^T θ‖₁。两种方式算出来就能直接对比边界误差。

完整复现代码如下,依赖 numpy、scipy、matplotlib:

import numpy as np from scipy.optimize import linprog from scipy.spatial import ConvexHull import matplotlib.pyplot as plt class Zonotope: def __init__(self, c, G): self.c = np.asarray(c, float).reshape(-1, 1) self.G = np.asarray(G, float) def storage_poly(P_max, E_max, E0, dt, T): """返回 np.array A, b 表示储能多面体约束 A x <= b""" A_p = np.vstack([np.eye(T), -np.eye(T)]) L = np.tril(np.ones((T, T))) * dt A_e = np.vstack([L, -L]) A = np.vstack([A_p, A_e]) b = np.r_[np.full(T, P_max), np.full(T, P_max), np.full(T, E_max - E0), np.full(T, E0)] return A, b def support_poly(A, b, theta): """多面体在 theta 方向上的支撑点,max theta^T x s.t. A x <= b""" res = linprog(-theta, A_ub=A, b_ub=b, bounds=[(None, None)] * len(theta)) if res.success: return res.x raise RuntimeError('支撑点求解失败') def support_zono(Z, theta): """奇诺多面体的支撑值 = theta^T c + sum |theta^T g_j|""" return float(theta @ Z.c.ravel() + np.abs(theta @ Z.G).sum()) def build_storage_zono(P_max, E_max, E0, dt, T): """简化单资源奇诺多面体:每时段功率生成器 + 能量耦合对角方向""" P_gen = np.hstack([np.eye(T) * P_max, -np.eye(T) * P_max]) e_dir = np.linspace(-0.5, 0.5, T)[:, None] * (E_max - E0) / (T * dt) G = np.hstack([P_gen, e_dir]) return Zonotope(np.zeros((T, 1)), G) def reduce_zonotope(Z, k=10, mode='inner'): U, s, _ = np.linalg.svd(Z.G, full_matrices=False) ratio = np.cumsum(s ** 2) / np.sum(s ** 2) k = min(k, len(s), int(np.searchsorted(ratio, 0.95) + 1)) U_k = U[:, :k] w = np.abs(U_k.T @ Z.G).max(axis=1) w *= 0.95 if mode == 'inner' else 1.05 return Zonotope(Z.c, U_k @ (w * np.eye(k))) # 场景参数 T, dt, N = 96, 0.25, 10 rng = np.random.default_rng(42) res_list = [storage_poly(0.5 + 0.5 * rng.random(), 2 + rng.random(), rng.random(), dt, T) for _ in range(N)] # 1) 真实聚合支撑:逐方向求各资源 LP 并求和 theta_list = [] for deg in np.linspace(0, 2 * np.pi, 72, endpoint=False): th = np.zeros(T) th[40] = np.cos(deg) th[41] = np.sin(deg) h_true = 0.0 for A, b in res_list: x = support_poly(A, b, th) h_true += th @ x theta_list.append((th, h_true)) # 2) 奇诺多面体聚合 + 缩减 zono_sum = Zonotope(np.zeros((T, 1)), np.zeros((T, 0))) for i, (P_max, E_max, E0) in enumerate([(0.5 + 0.5 * rng.random(), 2 + rng.random(), rng.random()) for _ in range(N)]): zono_sum = zono_sum.minkowski_sum(build_storage_zono(P_max, E_max, E0, dt, T)) zono_red = reduce_zonotope(zono_sum, k=8) # 3) 对比第 40、41 两个相邻时段的投影 h_red = [support_zono(zono_red, th) for th, _ in theta_list] h_true = [v for _, v in theta_list] rmse = np.sqrt(np.mean((np.array(h_red) - np.array(h_true)) ** 2)) print('支撑距离 RMSE =', round(rmse, 3), 'MW') print('缩减前生成器数 =', zono_sum.G.shape[1], '缩减后 =', zono_red.G.shape[1])

代码里最核心的是第 40、41 时段的支撑方向构造,这里用相邻时段做投影,能直接检验奇诺多面体对储能跨时段能量耦合的表达能力。运行时缩减函数会把 10 台设备拼接出的生成器从数百个压到 8 个,支撑距离的均方根误差通常在 0.05~0.15 MW 之间,取决于随机种子。

值得说明的是 build_storage_zono 是对单设备多面体的简化奇诺多面体近似,真实论文实现中这个环节会用约束奇诺多面体或迭代投影法替代,但聚合和缩减的主流程完全一致。先跑通这条链路,再替换资源级建模,是复现同类工作最稳妥的顺序。

4.2 结果中要重点看的三项指标

第一个指标是缩减前后支撑距离的 RMSE,建议按时段对画两条支撑曲线,能直观看到哪些角度方向误差偏大,通常出现在能量约束方向。第二个指标是生成器缩减前后数量比,工程上保留 5% 以内的生成器数是常见水平。第三个是二维投影面积,缩减后面积不应明显收缩(inner 模式收缩 5% 以内算正常),如果超过 20%,说明 k 太小或能量耦合方向被 SVD 截掉了,把 k 调大即可。

5. 聚合结果验证:包含率检查与高维投影调试技巧

5.1 用包含率做验收,不要只看支撑距离

支撑距离反映边界平均误差,但它不能回答最关键的问题:聚合结果里有多少指令是真实可执行的。对 inner 模式来说,需要统计从真实聚合域采样出的点中,有多大比例落在缩减奇诺多面体内部。代码里用一个线性规划可行性检查逐点判断:

from scipy.optimize import linprog def point_inside_projected_zono(q, Z, R): """判断投影点 q 是否在 R*Z 内:RG z = q - R c, -1 <= z <= 1""" RG = R @ Z.G rhs = (q.reshape(-1, 1) - R @ Z.c).ravel() A_eq = RG n = RG.shape[1] bounds = [(-1, 1)] * n res = linprog(np.zeros(n), A_eq=A_eq, b_eq=rhs, bounds=bounds) return res.success # 从真实投影多边形的凸包中随机采样(此处示意三角形采样,完整版用 Delaunay) # 省略采样实现,核心是每个点调用 point_inside_projected_zono 并统计成功率

包含率低于 95% 时优先怀疑两个位置:一是缩减函数里的 0.95 缩放系数对边界场景不够保守;二是 SVD 截断把能量方向削掉了,解决方法是把 k 提高到 12~15 或改用带约束的线性规划求生成器长度。

5.2 高维投影的二维调试法

直接看 96 维空间里的聚合结果毫无意义,我一般固定两个相邻时段构造投影矩阵 R,把高维集合画成二维多边形,再叠加真值支撑点。调试时依次试三个投影:相邻时段、间隔 12 时段、间隔 48 时段。相邻时段投影反映功率爬坡与能量连续性,间隔 48 时段反映跨日能量转移能力。如果某个投影方向误差特别大,对应时间尺度上的耦合关系被压缩得过度,就应该加大该时段的权重,或在生成器缩减前对 G 按列做加权。这个方法比单纯调 k 或缩放系数更快定位问题。

奇诺多面体的可行域聚合在工程上的价值是牺牲少量精度换取计算可行性,验证工作的核心就是把这个“少量精度”量化清楚。支撑函数和包含率两个工具配合,能把近似的边界误差和可执行风险都摆在明面上,后续无论接入市场申报还是实时调度,都按这套校验结果决定是否启用。

本文还有配套的精品资源,点击获取

需要专业的网站建设服务?

联系我们获取免费的网站建设咨询和方案报价,让我们帮助您实现业务目标

立即咨询