企业官网建设流程全解析

我最早接触含冰蓄冷的冷电联供微网，是在一个制冷季的能源审计项目里。白天办公楼冷负荷爆表，正好撞上电价高峰段；夜里空调没需求，电价反而跌到谷段的三分之一。当时脑子里蹦出来的第一个念头就是：要是能在夜里把电变成冰存起来，白天再用这些冰扛冷负荷，电费不就省下来了？等真正把它写进微网经济优化模型里才发现，冰蓄冷远不止“搬移冷量”这么简单——它会牵动燃气轮机出力、电网购电、电制冷机投停、蓄冰槽调度的一连串连锁反应。这篇文章就是把我在MATLAB里完整跑通“含冰蓄冷装置的冷电联供型微网经济优化运行”的过程、模型思路和踩过的坑整理出来，适合正在做综合能源、微电网优化方向的研究生，以及做园区能源规划的工程师参考。

这类课题的核心其实可以拆成三个问题：系统里有哪些设备、冷和电怎么流动；经济优化的目标函数和约束条件怎么用数学语言表达；以及如何用MATLAB把模型建起来、求解出来、把结果讲清楚。下面我就按这个逻辑一步步展开，最后把代码实现和典型日算例结果一起贴出来，保证你照着做能跑出图。

1. 先把这个系统的能量流彻底捋清楚

1.1 微网里都有谁，谁在供电谁在供冷

含冰蓄冷的冷电联供微网，典型的设备清单大概是这么一套：燃气轮机（或者燃气内燃机）、余热制冷机组（通常是溴化锂吸收式）、电压缩制冷机、冰蓄冷装置，再加上电负荷和冷负荷这两个“必须满足的口子”。有些配置里还会加燃气锅炉作为备用热源，但这个课题如果只关心冷电联供，燃气锅炉可以先不放进模型。

燃气轮机烧天然气发电，电优先供给微网内的电负荷；发完电剩下的烟气余热不能浪费，送进溴化锂机组去制冷，这部分冷量就是“联供”的核心产出。电制冷机的角色比较灵活：白天冷负荷高的时候它可以直接制冷补缺口，夜里电价便宜的时候它又可以把冷量以冰的形式存在蓄冰槽里。而冰蓄冷装置本身不产冷，它只是把“夜间便宜的电”转化成“白天能用的冷”，相当于在时间轴上做能量的搬运。

把这套系统画成能量流来看：燃气轮机的出口有两条支路，一条是电母线，一条是烟气余热母线；电母线上挂着电负荷、电制冷机，以及一个可选的“向电网购电/售电”接口；冷母线上挂着溴化锂机组、电制冷机、融冰出口和冷负荷。优化调度做的事情，就是每个调度时段决定燃气轮机发多少电、电网买多少电、电制冷机产多少冷且多少直接供冷多少去蓄冰、蓄冰槽融多少冰出来，最终让电和冷两条母线都刚好平衡。

这里我建议初学者先把“冷母线”和“电母线”分开画，两个口子的平衡约束分开写。很多第一次做这个课题的人，容易把电制冷机的耗电功率和它的产冷量混在一起，一写约束就乱套。记住一条：电平衡只算功率，冷平衡只算冷量，电制冷机就是它们之间的“转换器”，用COP把两者勾起来。

1.2 冰蓄冷为什么能省钱，它的“天时”在哪

冰蓄冷省钱的根本逻辑，不是设备效率，而是时间套利。理解这件事特别重要，因为它决定了你后面积分时电价、冷负荷曲线时该怎么思考。

举个生活的例子：假设你所在的地区白天电价1.2元/度，夜间谷段电价0.4元/度，白天开空调制冷需要的电费是晚上的三倍。冰蓄冷相当于你在晚上用电把水冻成冰，白天冰块融化的冷量来吹空调，等于把白天的用电量“挪”到了晚上交费。表面上总耗电量没怎么变，但电费结构变了，这就是它经济性的第一来源。

在含燃气轮机的冷电联供微网里，经济性会比“单纯峰谷套利”复杂一层。因为燃气轮机烧燃气发电，也要算燃料成本；它余热制冷还能顶掉一部分电制冷机的耗电。所以系统最优调度不是简单地“谷段蓄冰、峰段融冰”，而是要全局打量：某个时段到底是买电网的电更划算，还是自己烧气发电更划算；发出来的电是给电负荷用还是给电制冷机用；冷负荷缺口是用溴化锂余热供、还是电制冷直供、还是把蓄冰槽里的冰放出来。这些都是优化模型需要解决的问题。

还有一点容易被忽略：并不是所有冷负荷都适合用冰蓄冷来扛。冰蓄冷装置有一个蓄冰容量上限，蓄满了就是蓄满了，放空了就是放空了。而且融冰出力的速率也有限制，不能指望一个小时内把整夜的冰全部融完。所以它更像一个“蓄水池”，有容积约束、有进出水速率约束，调度模型必须把它当储能设备来建模，而不是一个可以任意使用的“免费冷源”。

2. 经济优化模型究竟在优化什么

2.1 目标函数：一天下来究竟花多少钱

经济优化运行的“经济”二字，必须落到一个可量化的函数上。最常见的做法是取一个典型日，把全天分成T个调度时段（通常取24个时段，每小时一个点），目标函数是这24个小时里微网的总运行成本最小化。

总运行成本一般包括三块：从电网购电的费用、燃气轮机消耗天然气的燃料费、以及各设备的运维费用。写出来大概是这样的形式：

min ∑(购电电价 × 购电功率 + 气价 × 燃气耗量 + 运维系数 × 出力)

其中购电电价是分段函数，峰、平、谷各时段电价不一样，这在MATLAB里可以用一个长度24的向量直接表示，每一个元素就是对应时段的电价。燃气耗量跟燃气轮机出力之间的关系，可以用一个线性的热耗率表达式来近似：燃气耗量 = 燃气轮机出力 / 发电效率，再换算成燃料费。运维费用通常按设备的运行出力乘以一个很小的单位运维成本系数来估算，比如电制冷机每产出1kW冷量收0.01元、燃气轮机每发1kW电收0.02元，这部分费用虽然数额不大，但加入目标函数之后能让优化结果更贴近实际，避免出现设备出力频繁跳变的毛刺。

有些模型还会考虑微网向电网售电的收入，也就是如果燃气轮机发多了电、用不完，可以卖给电网赚钱。这种情况下目标函数里就要减掉卖电收入，但大多数课题为了聚焦冰蓄冷的调度逻辑，会先假设微网不从电网卖电，或者只允许从电网买电不允许反向送电。建议你入门阶段先做成“只买不卖”，把主框架跑通了再放开。

2.2 约束条件：电力平衡、冷量平衡和蓄冰槽的“脾气”

目标函数定了之后，约束条件就是模型的骨架。缺了约束的优化就是一个无约束最大化利润问题，很荒谬。对这套系统来说，至少有以下几类约束必须写清楚。

第一是电功率平衡约束。每个时段，燃气轮机出力加上从电网购电的功率，必须等于电负荷加上所有用电设备的耗电功率。用电设备包括电制冷机（注意它既可能直接制冷也可能蓄冰，两部分耗电都要计入），以及其他辅助设备。这条约束是等式约束，物理含义是“微网内部电能不能凭空产生或消失”。

第二是冷功率平衡约束。每个时段，溴化锂机组产冷量、电制冷机直供冷量、蓄冰槽融冰供冷量三者之和，必须等于冷负荷。同样也是等式约束。这里溴化锂机组的产冷量跟燃气轮机的余热息息相关，一般用热电比或者余热回收效率把它和燃气轮机出力耦合成一个线性关系。

第三是蓄冰槽的动态约束。这是整个模型中最有“储能味”的一条。蓄冰槽的蓄冷状态可以用SOC（State of Charge，荷冷状态）来表示，它像一个水桶：本时段结束时桶里的水量，等于上一时段剩下的水量，加上本时段蓄进去的冷量（蓄冰），减去本时段放出来的冷量（融冰）。写成离散形式就是：

SOC(t+1) = SOC(t) + (蓄冰功率 - 融冰功率) × 时段时长 / 蓄冰槽容量

SOC不能小于0也不能大于1，这对应蓄冰槽不能“透支”也不能“溢出”。另外，蓄冰槽的蓄冰速率和融冰速率都各自有上限，取决于制冷机和融冰换热器的实际能力。对于完整的日调度，通常还会要求SOC在一天开始和一天结束时保持一致，代表这个调度策略可以循环运行。

第四是设备出力上下限约束和爬坡约束。燃气轮机有最小稳定出力和最大额定出力，电制冷机的制冷能力也有上下限。爬坡约束描述的是设备从一个时段的出力跳到下一个时段时，变化速度不能太快，比如燃气轮机每15分钟最多只能爬升百分之多少，这对燃机实际运行很重要。

第五是互斥约束。同一时刻，蓄冰槽在蓄冰和在融冰是两种矛盾状态，理论上不应该同时发生。这需要引入二进制变量，比如用0-1变量y表示蓄冰工况：y=1时只能蓄冰、不能融冰，y=0时只能融冰、不能蓄冰。加了二进制变量，整个规划问题就变成了混合整数线性规划MILP，求解规模会变大，但这也是这类课题的标准做法。

2.3 关键参数的计算与单位陷阱

参数处理是实操中第一个“翻车高发区”。我见过太多同学把冷负荷单位搞混，结果优化出来的结果全是“天文数字”。这里必须提一下几个最常见的单位坑。

冷量单位有kW和RT（冷吨）两种，二者之间的换算是1RT约等于3.517kW。很多文献里冷负荷曲线用的是RT，电负荷用的是kW，如果直接放进同一个模型里算，冷平衡约束就废了。我习惯的做法是在数据读入MATLAB之后，统一把冷负荷单位转换成kW再参与计算，所有设备参数也都转成kW制，最后出图时再根据需要标回RT。

蓄冰槽的容量常见有两种标法：一种是用蓄冷量表示，比如“××RTh”或“××kWh”，另一种是用蓄冰量表示，比如“××吨冰”。如果是后者，换算更麻烦一点，因为冰的融解热是334kJ/kg，1吨冰全部融化大约可以吸收93kWh左右的热量。实际建模时，我建议直接以蓄冰槽的冷量容量为基准定义SOC，这样最简单。

COP的取值也要分场景。同样是那台电制冷机，直接制冷的COP和制冰工况的COP通常不一样，制冰工况的冷凝温度更低、效率大约会打个八折。模型里最好拆成两个参数：COP_direct和COP_ice，蓄冰功率除以COP_ice才是它的耗电功率。很多简化模型干脆只用一个COP，误差其实不小，尤其当蓄冰量占比较大时，算出来的夜间耗电量会偏低，省电费效果会被高估。

还有一个容易忽略的是分时电价边界时段的处理。比如峰段从8:00开始，那么7:59和8:00这两个时段如果分别落在谷段和峰段，电价会出现一个突跳。优化模型会强烈倾向于把大量蓄冰安排在谷段末尾，这是合理的；但如果蓄冰速率上限设置得过大，模型可能安排极端出力，出现实际设备无法响应的调度指令。所以蓄冰速率上限要按设备的真实能力来约束，不能拍脑袋给一个很大的数。

3. MATLAB求解路线怎么选

3.1 方案一：YALMIP + Cplex 精确求解

做混合整数线性规划，我在MATLAB里最常用的组合是YALMIP加外部求解器。YALMIP是一个建模工具箱，它最大的价值是让你把优化问题写成接近数学公式的代码，而不用关心底层求解器的具体调用方式。Cplex或者Gurobi扮演的是求解引擎的角色，负责真正把MILP解出来。两者配合起来，写代码的效率远高于手写linprog矩阵。

为什么推荐这条路？因为这个课题的模型里天然有二进制变量，用YALMIP可以非常自然地写出“蓄冰和融冰互斥”这类逻辑约束。如果不用YALMIP，你自己用ma tlab的intlinprog来建模，就必须把所有的约束条件手动整理成矩阵形式Ax≤b，而且要把二进制变量和连续变量拼成一个长向量，变量一多，下标映射能把人绕晕。YALMIP的典型写法是：

T = 24; x = sdpvar(T, 1); % 连续变量 y = binvar(T, 1); % 二进制变量 Constraints = []; for t = 1:T Constraints = [Constraints, x(t) >= 0, x(t) <= 100]; end Cost = sum(x); optimize(Constraints, Cost, sdpsettings('solver', 'cplex'));

跑完用value()函数取出各变量的数值，然后直接plot，非常顺手。需要提醒的是，Cplex和Gurobi都是商业求解器，但如果是在校学生，可以通过学校许可证或者求解器官方申请学术版授权。实在装不上，也可以用YALMIP内置支持的免费求解器如CBC、GLPK。免费的求解器速度会差一些，但这个课题的规模不大，一般也能在可接受的时间内跑完。

3.2 方案二：启发式算法兜底

有些朋友的环境里装不了外部求解器，或者导师坚持要求用智能优化算法，这时会考虑粒子群、遗传算法这类启发式方法。我不拦着，但必须先说清楚它的代价。

启发式算法最麻烦的是处理约束。MILP里那些等式约束（电平衡、冷平衡）和不等式约束，在粒子群算法里不能像YALMIP那样直接写进去。常见的做法是把约束转成罚函数加到目标函数里，或者设计满足约束的编码方式。罚函数法说白了就是“如果违反约束，就在目标函数里加一个很大的惩罚值”，但惩罚系数怎么取是个玄学：取小了，结果违反约束；取大了，惩罚项直接盖过真实成本，优化过程变成在瞎找。我见过不少代码，罚函数系数拍脑袋设了个10000，结果目标函数全是罚出来的数，完全失去经济意义。

另一个问题是重复性差。启发式算法有随机性，每次跑的最终结果都不一样，除非设定随机数种子。写论文做对比实验时，“优化结果不稳定”这个缺点非常恼人。所以我的建议是：如果条件允许，优先用精确求解器跑MILP；启发式算法可以作为对比方法，在论文里验证智能算法的有效性，但不建议一开始就陷进去。

3.3 核心代码结构与关键片段

下面这段是我整理的一套可运行的YALMIP关键代码骨架，覆盖了变量定义、约束书写和求解调用。参数部分我统一写成一个结构体，方便后续修改。

%% 数据准备（假设每个变量都是长度为T的列向量） T = 24; dt = 1; % 调度时段数和步长 P_load = [...]'; % 24h电负荷，kW Q_load = [...]'; % 24h冷负荷，kW（注意单位统一） c_buy = [谷段电价; 峰段电价; ...]'; % 24h购电价 % 设备参数 cap_ice = 500; % 蓄冰槽容量，kWh rate_charge = 150; % 蓄冰速率上限，kW rate_melt = 200; % 融冰速率上限，kW COP_ec = 3.5; % 电制冷机直供COP COP_ice = 2.8; % 电制冷机制冰工况COP COP_br = 1.2; % 溴化锂余热制冷系数 eta_gt = 0.35; % 燃气轮机发电效率 gas_price = 2.6; % 天然气价，元/m3 P_gt_max = 500; P_gt_min = 50; % 燃气轮机出力上下限 Q_ec_max = 600; % 电制冷机最大产冷量 %% 变量定义 P_gt = sdpvar(T,1); % 燃气轮机出力 P_buy = sdpvar(T,1); % 购电功率 Q_ec_direct = sdpvar(T,1); % 电制冷机直供冷量 Q_ec_charge = sdpvar(T,1); % 电制冷机用于蓄冰的冷量 Q_melt = sdpvar(T,1); % 融冰供冷量 soc = sdpvar(T+1,1); % 蓄冰槽SOC y = binvar(T,1); % 蓄冰/融冰互斥标志 %% 约束 Constraints = []; for t = 1:T % 电平衡 Constraints = [Constraints, P_gt(t) + P_buy(t) == P_load(t) + ... Q_ec_direct(t)/COP_ec + Q_ec_charge(t)/COP_ice]; % 冷平衡 Constraints = [Constraints, COP_br*P_gt(t) + Q_ec_direct(t) + Q_melt(t) == Q_load(t)]; % 设备上下限 Constraints = [Constraints, P_gt_min <= P_gt(t) <= P_gt_max]; Constraints = [Constraints, 0 <= P_buy(t) <= 1000]; Constraints = [Constraints, 0 <= Q_ec_direct(t) <= Q_ec_max]; Constraints = [Constraints, 0 <= Q_ec_charge(t) <= rate_charge]; Constraints = [Constraints, 0 <= Q_melt(t) <= rate_melt]; % 互斥约束 Constraints = [Constraints, Q_ec_charge(t) <= rate_charge*y(t)]; Constraints = [Constraints, Q_melt(t) <= rate_melt*(1-y(t))]; % 蓄冰槽动态 Constraints = [Constraints, soc(t+1) == soc(t) + ... (Q_ec_charge(t) - Q_melt(t))*dt/cap_ice]; Constraints = [Constraints, 0 <= soc(t+1) <= 1]; end % 循环周期约束 Constraints = [Constraints, soc(1) == soc(T+1)]; %% 目标函数与求解 Cost = sum(c_buy.*P_buy*dt + ... gas_price*P_gt*dt/eta_gt/10 + ... % 天然气耗量换算，注意气热值单位 0.02*P_gt + 0.01*(Q_ec_direct+Q_ec_charge)/COP_ec*dt); ops = sdpsettings('solver','cplex','verbose',1); result = optimize(Constraints, Cost, ops); if result.problem == 0 P_gt_opt = value(P_gt); P_buy_opt = value(P_buy); soc_opt = value(soc); else disp('求解失败，检查约束和求解器设置'); end

这段代码我尽量保持了精简，但核心要素都齐全。运行起来之后，你可以把P_gt_opt、soc_opt这些结果画成柱状图或阶梯图，观察各设备的出力时序。如果遇到具体报错，后面第5节我会分析常见问题。

4. 典型日算例结果怎么看

4.1 场景与参数设定

为了让你对结果有个直觉，我说一个自己调过的典型算例。假设一个商业园区，冷负荷高峰出现在下午14:00到17:00，峰值大约是800kW；电负荷跟冷负荷走势相近但耦合不强，峰值约600kW。分时电价按常见的三段式：谷段0:00-8:00电价为0.4元/kWh，平段8:00-11:00和17:00-21:00为0.8元/kWh，峰段11:00-17:00和21:00-24:00为1.2元/kWh。燃气轮机最大出力500kW，最小稳定出力50kW，发电效率35%；溴化锂余热制冷系数取1.2；电制冷机最大产冷量600kW，COP取3.5，制冰工况COP取2.8；蓄冰槽容量500kWh，蓄冰速率上限150kW，融冰速率上限200kW。天然气价按2.6元/m³，再按天然气热值折算单位发电成本。

这个场景我跑出来的典型优化结果有几个鲜明特征。第一，蓄冰槽的SOC曲线呈现出非常典型的“夜间蓄、白天放”形态：从凌晨1点左右开始SOC一路爬升，到早上8点前达到接近1，然后从上午10点左右开始持续下降，到傍晚冷负荷高峰期基本放空。第二，燃气轮机并非全天满发，而是在电价较高的峰段加大出力，在谷段则压低出力甚至接近最小出力，因为谷段买电网电更划算。第三，电制冷机在谷段大量运行，一部分冷量直供夜间的基数冷负荷，大部分冷量用于制冰蓄存；到了白天峰段，电制冷机的直供冷量会被压缩，冷负荷主要靠融冰和余热制冷顶上。

4.2 优化结果怎么解读，对比数据才有说服力

光看一个优化结果没啥说服力，做研究必须要有对比。最直接也最常用的对比方案，是设一个“无冰蓄冷”的对照组：把冰蓄冷相关变量全部置0，蓄冰槽容量设为0，其他参数保持不变，重新求解一个不含冰蓄冷的冷电联供微网优化模型。两个模型跑出来的总运行成本一对比，冰蓄冷的经济价值就出来了。

在我这个算例里，含冰蓄冷方案的全天总运行成本大约是4800元，无冰蓄冷方案大约是5500元左右，节约比例在10%到15%。这个数字会让你直观感受到“移峰填谷”的能量有多大。你还可以继续做敏感性分析：拉大峰谷电价差，看节约比例怎么变；或者把蓄冰槽容量从300改到800，看经济性提升是否边际递减。这类分析放在论文里就是非常扎实的章节。

结果图的呈现上，我建议至少画三张图。第一张是蓄冰槽SOC曲线，一眼看出蓄放规律；第二张是冷负荷供需堆叠图，把溴化锂供冷、电制冷直供、融冰供冷三个部分堆叠起来，跟冷负荷曲线对整齐，谁在什么时段出力一目了然；第三张是电功率平衡图，把燃气轮机出力、电网购电、电负荷和电制冷耗电画在一起。这几张图画完之后，整个系统的调度策略就非常立体了。

5. 实操中常见的问题与排查技巧

5.1 求解器与MATLAB代码层面的报错

我在调试这套模型时遇到过不少报错，挑几个高频的写在这里供你排查。

“No suitable solver found”是YALMIP最常见的报错。它表示求解器没有正确配置到YALMIP的路径里。解决方法是检查Cplex或Gurobi的MATLAB接口是否安装，然后用yalmiptest命令测试一下哪些求解器可用。对于MILP问题，YALMIP需要有支持整数规划的求解器；如果你只有linprog，那确实解不了带binvar的模型。一个应急办法是把binvar去掉，假设蓄冰和融冰可以同时发生，模型退化成LP，用自带linprog也能跑，但结果在物理上不严格。

“Index exceeds array bounds”这种报错大多是变量维度不匹配。YALMIP里如果soc定义成(T+1,1)，而循环里误写成soc(t)，就很容易越界。另外，约束中两个变量相加时，要注意它们都是列向量还是行向量。MATLAB对1×24和24×1的隐式扩展有时能算，但结果会莫名其妙。我的习惯是所有时间序列变量统一定义为T×1列向量，绝不混用。

“Infeasible problem”提示无可行解是最让人头疼的。我第一次跑就遇上了，后来发现是SOC的周期约束和初始SOC设置矛盾了：我既设了soc(1)=0.5，又要求soc(T+1)=soc(1)，这本来没问题，问题是我在凌晨安排了很大的蓄冰量，导致soc(1)=0.5的约束限制了初值偏高，无法在24小时内完成一次完整的蓄放循环。解决方法是先跑一次不带SOC周期约束的模型，观察末尾SOC自然落点，再来决定周期约束怎么写，或者干脆让soc(1)=soc(T+1)=0，表示一天从空槽开始也相对合理。

5.2 模型与工程层面的坑

单位不一致是最隐蔽的坑。冷负荷曲线如果是从RT标定的，而蓄冰槽容量是用kWh标定的，你必须在数据准备阶段就统一成kW和kWh。我建议在代码开头加一个单位转换参数，比如RT_TO_KW = 3.517，所有冷量数据乘完再进模型，同时用注释写明转换关系，免得过几天自己都忘了。

蓄冰槽的自损耗也不要完全忽略。实际蓄冰槽一天下来会有5%到10%的冷量损失，如果你完全不建模，优化结果会过于乐观。最简单的方法是在SOC动态约束里加一个衰减系数，比如SOC(t+1) = 0.98×SOC(t) + 入流 - 出流，虽然只是粗略近似，但结果的可信度提升不少。

电制冷机同时承担直供和制冰时，COP分开取值的细节前面已经强调过。这里再补充一个实操技巧：如果你发现优化结果中电制冷机蓄冰速率频繁触及上限，说明蓄冰槽容量相对冷负荷峰值偏大，经济上出现了“为了蓄冰而蓄冰”的边际效益递减。此时可以做一个蓄冰槽容量敏感性扫描，找出边际收益骤降的拐点，这个拐点就是合理的容量配置。

最后提醒一点，关于爬坡约束。如果燃气轮机的爬坡约束设置得太紧，求解时间会明显变长，因为二进制变量和连续变量之间的耦合更紧了。我一般先不加入爬坡约束，快速验证整个模型框架无误后，再把它加回去做完整仿真。这样能大幅缩短前期调试时间。

我个人做完这个课题后最大的体会是：冰蓄冷装置的引入，真正改变了微网调度的“相位”——把电力需求从白天挪到黑夜，把冷量供给从即时生产变成了“先储后放”。MATLAB在这个课题里最大的价值，是借助YALMIP把模型和算法解耦，让研究者的精力集中在物理建模和结果分析上，而不是陷入求解器的矩阵地狱。最后再分享一个小技巧：拿到任何冷电联供微网优化模型，先跑一个不蓄冰的纯电制冷场景，把电平衡、冷平衡和负荷曲线彻底验证对了，再一步步加入余热制冷、冰蓄冷、互斥约束。按这个顺序递进，你会发现自己调试模型的效率至少翻一倍。