这两年做综合能源系统方向的人越来越多了。不管是高校课题组做园区级能源规划,还是设计院做区域供能方案,基本都绕不开一个基础问题:把电、热、气这几个系统放在一起算潮流。纯电网的潮流计算大家都会,热网水力计算也好办,但是一旦牵扯到多能耦合,电网和热网之间互相传递功率,电气热能流就变成了一个非线性强耦合问题。身边不少人卡在这一步,程序要么不收敛,要么结果明显不合理,最后只能回头去翻论文、抄代码。
这篇博文就是围绕项目标题“计及多能耦合的区域综合能源系统电气热能流计算研究(Matlab代码实现)”来展开。我会把整个研究的出发点、数学模型、Matlab实现思路、算例结果和排查经验都讲一遍,重点是告诉大家实际跑代码时会遇到什么问题,哪些地方教材上不会写。无论你是刚接触综合能源的研一学生,还是想快速搭建区域能源系统仿真平台的工程师,这篇文章应该都能给你省不少时间。
1. 项目整体设计与需求拆解
1.1 为什么综合能源系统必须做多能耦合潮流计算
传统电力系统潮流计算只关心电功率平衡,节点上无非是发电机、负荷和线路。综合能源系统不一样,电、热、气三套网络在物理上交织在一起:燃气轮机烧天然气发电,余热进热网;电锅炉吃电产热,直接改变电网负荷分布;热泵从电网取电,向热网注入热量。这三个网络谁也不独立,必须联合求解才能得到稳态运行点。
如果忽略耦合关系,分别算电网和热网,会出大问题。比如电锅炉在电网侧是个大负荷,在热网侧是个热源,你单独算电网时把它当恒定负荷,单独算热网时把它当恒定热源,两边都假设对方不变,可实际上电网电压变化会影响电锅炉的实际电功率,电锅炉的电功率又直接决定热出力,这种强反馈关系只有联立求解才能如实反映。这也是这个项目标题里“计及多能耦合”六个字的含义所在。
1.2 Matlab实现方案的整体选型考量
Matlab在这个领域依然是首选环境。综合能源系统潮流计算本身需要处理矩阵运算、稀疏存储和非线性迭代,Matlab的矩阵操作和内置的求解器让原型搭建效率很高,比C++和Python都更省事。尤其是做学术研究,很多时候需要反复修改网络拓扑和参数,Matlab的脚本式开发方式改起来很快。
代码实现上我采用了模块化分层结构,每个网络独立成函数,耦合组件的模型单独封装,外层用交替迭代法做全局收敛。这样设计的好处是可以任意替换某个子网的求解器,比如把电网牛拉法换成改进的快速分解法,或者把热网节点法换成面向环路的水力求解器,都不影响其他模块。
说句题外话,工程上有些大软件包像EBSILON和TRNSYS也能做多能系统仿真,但它们的潮流计算往往内置在黑盒里,用户很难干预算法细节,也很难跟自己的数据接口对齐。Matlab代码是自己可控的,模型透明,适合做方法研究。
1.3 测试算例的方案设计
考虑到通用性,我把测试系统设计成一个典型区域综合能源系统:电网用IEEE 33节点配电网做基础,热网用一个6节点辐射状管网,气网用6节点输气网络,通过CHP机组、电锅炉和燃气锅炉三个耦合设备连接起来。这个规模不大不小,既能体现多能耦合的完整机理,又不至于让代码运行时间太长。
算例设计上我刻意留了一些“坑”:热网负荷与电网负荷峰谷不完全同步,CHP电出力上限约束在某个值附近,气网管道压降比较大,天然气供应处于紧张状态。这些设计在后面的结果分析中都会体现出来,能比较全面地考察多能耦合计算方法的收敛性和物理合理性。
2. 电气热多网数学模型与核心原理
2.1 配电网潮流计算的数学模型
电网侧采用牛顿-拉夫逊法,这是经典成熟的做法。对每个PQ节点,功率不平衡量方程为:
ΔPi = Psp,i - Ui∑Uj(Gij cosθij + Bij sinθij) = 0
ΔQi = Qsp,i - Ui∑Uj(Gij sinθij - Bij cosθij) = 0
其中Psp和Qsp是节点给定的有功和无功注入,Gij和Bij是节点导纳矩阵的实部虚部。求解时先建立节点导纳矩阵Ybus,然后反复迭代修正电压幅值和相角,直到功率不平衡量小于阈值。
这个模型本身没什么新鲜的,但接入多能耦合设备后有两个细节需要注意。第一,CHP机组在电网中通常作为PQ节点或PV节点处理,取决于控制策略,如果是“以热定电”,CHP的电出力由热负荷决定,此时应作为PQ节点;如果是“以电定热”,则电出力固定,热出力跟随变化。第二,电锅炉在电网侧就是纯负荷,但负荷值不是恒定的,因为它取决于热网侧的热负荷需求,需要在耦合迭代中实时更新。
2.2 热力管网水力与热力计算的数学建模
热网的计算比电网麻烦得多,因为它包含水力计算和热力计算两层。水力计算的核心是流量连续性方程和环路压降方程:
对任意节点,流入流量等于流出流量加节点流量消耗;
对任意管道回路,所有管道压降之和等于零。
管道压降通常用 Darcy-Weisbach 公式或经验公式计算,涉及摩擦系数、管径、长度和流量。传热计算则是基于热功率平衡,对每个热负荷节点有:
φ = cp·m·(Ts - Tr)
其中φ是热负荷功率,cp是水的比热容,m是流经节点的质量流量,Ts和Tr分别是供水和回水温度。
实际算热网时,我推荐采用图论中的关联矩阵方法建方程。节点-支路关联矩阵A描述拓扑关系,环路-支路关联矩阵B描述环网关系,联立求解各管道流量,再以节点温度递推的方式求解管网温度分布。这样做代码结构清晰,扩展环路或增加节点都只是改矩阵的事,不需要重写方程。
2.3 天然气网络与耦合设备模型
天然气网的水力计算使用Weymouth稳态方程,描述管道流量与两端压力平方差的关系。这个方程也是非线性的,求解方式和电网潮流类似,只是状态变量是节点压力而不是电压,需要单独写一套牛顿迭代。
耦合设备是整套计算中的关键纽带。我建立了三种设备模型,第一是CHP机组,它有电出力P_e、热出力Q_h和燃料耗量F三者的映射关系,我用简化的线性模型:
Q_h = C_m·P_e
F = P_e / η_e
其中C_m是热电比,η_e是发电效率。这种简化在区域供能场景下足够用,但如果研究变工况还要替换成更精细的曲线模型。
第二是电锅炉,效率取0.95,热出力等于电功率乘效率。第三是燃气锅炉,直接消耗天然气产热,效率按0.9处理。这三种设备的模型都封装成独立函数,输入是需求侧参数,输出是耦合变量,外部循环直接调用即可。
3. Matlab代码实现与迭代求解流程
3.1 代码整体架构与数据结构设计
代码采用两个层次的迭代结构。内层是各个子系统的独立潮流求解,外层是全局多能耦合收敛迭代。我设计了七个核心函数:
buildGridYbus.m:构建配电网节点导纳矩阵solvePowerFlowNR.m:电网牛顿-拉夫逊法潮流求解buildHeatTopology.m:构建热网拓扑关联矩阵solveHeatNetwork.m:热网水力-热力联合求解solveGasNetwork.m:天然气网稳态潮流求解coupleDeviceModel.m:耦合设备多能转换模型main_multienergy.m:主程序,统筹全局迭代
数据结构上我统一使用结构体变量。比如电网节点参数放在gridNode结构中,包含编号、类型、有功负荷、无功负荷、电压初值等字段;热网管道数据放在heatPipe结构中,包含首末节点号、管长、管径、粗糙度等信息。用结构体的好处是传参方便,调试时一眼就能看出哪个字段缺失,比散落的数组变量强得多。
3.2 交替迭代法与统一求解法的路径选择
多能耦合潮流求解有两条技术路线:一是把电网、热网、气网方程统一组装成一个大规模非线性方程组,用整体雅可比矩阵求解,即统一求解法;二是各网络独立求解,通过耦合变量在外层来回传递,直至收敛,即交替迭代法。
统一求解法理论收敛性好,但工程实现难度大。因为电网潮流状态量是电压幅值和相角,热网状态量是节点压力和温度,气网状态量是节点压力,物理量纲不同、数值尺度差异巨大,直接拼装会导致雅可比矩阵条件数非常差,极难收敛。我试过一次,全场几乎没有正常收敛的例子。
所以这个项目最终选了交替迭代法。外层循环的流程是:先假定耦合设备的功率初值,将CHP电出力和电锅炉电功率作为电网节点注入,求解电网潮流;再把CHP热出力和电锅炉热出力作为热网热源,求解热网潮流;接着把CHP和燃气锅炉的耗气量作为气网节点流量,求解气网潮流;最后用新得到的网络状态更新耦合设备功率,检查前后两轮功率差是否满足收敛条件。
交替迭代法虽然收敛速度不如统一法,但胜在稳定可靠,而且每个子网都能应用各自领域成熟的求解技术,代码复用性好。
3.3 电网潮流求解核心代码实现
电网潮流用牛拉法实现。核心代码如下:
function [V, theta, iter] = solvePowerFlowNR(Ybus, Ssp, V0, tol, maxIter) % 输入: Ybus 节点导纳矩阵, Ssp 复功率注入, V0 电压初值 % 输出: V 节点电压幅值, theta 节点相角, iter 迭代次数 nb = length(V0); V = abs(V0); theta = angle(V0); for iter = 1:maxIter % 计算节点功率注入 Vc = V .* exp(1j*theta); Icalc = Ybus * Vc; Scal = Vc .* conj(Icalc); % 功率不平衡量(平衡节点不参与) dP = real(Ssp(2:end) - Scal(2:end)); dQ = imag(Ssp(2:end) - Scal(2:end)); dPQ = [dP; dQ]; if max(abs(dPQ)) < tol break; end % 构造雅可比矩阵 J = calJacobian(Ybus, V, theta); % 求解修正方程 dx = -J \ dPQ; % 更新状态量 theta(2:end) = theta(2:end) + dx(1:nb-1); V(2:end) = V(2:end) + dx(nb:end); end end这里有一个很关键的细节:构造雅可比矩阵时,一定要区分对角块和非对角块的公式。对角块上有注入功率对节点电压的偏导,还包含该节点自身导纳分量,很容易写错。我建议参考标准电力系统分析教材上的公式逐个元素推导一遍,不要直接用数值差分代替雅可比,那样不仅速度慢,而且在电压无功强烈的配电网中容易引入数值误差。
3.4 热网联合求解核心代码实现
热网求解是这套代码里最容易出错的部分。我的实现思路是先做水力计算,求解各管道流量,再做热力计算,求解节点供回水温度。
function [Ts, Tr, m_flow] = solveHeatNetwork(heatNode, heatPipe, Qin, Toutdoor) % 输入: heatNode 节点参数结构, heatPipe 管道参数结构 % Qin 各热源注入热功率, Toutdoor 室外温度 % 输出: Ts 供水温度, Tr 回水温度, m_flow 各管道质量流量 nn = length(heatNode); np = length(heatPipe); % 步骤1: 建立节点-支路关联矩阵 A(nn x np) A = zeros(nn, np); for k = 1:np i = heatPipe(k).from; j = heatPipe(k).to; A(i,k) = 1; A(j,k) = -1; end % 步骤2: 求解节点净流出流量平衡 % 根据热负荷需求和设计供回水温差,计算节点流量需求 dT_design = 20; % 设计温差 20度 cp = 4187; % 水的比热容 J/(kg·K) for i = 1:nn if heatNode(i).load > 0 heatNode(i).m_req = heatNode(i).load * 1000 / (cp * dT_design); else heatNode(i).m_req = 0; end end % 步骤3: 简化辐射状管网流量分配(从热源节点逐级向外推) m_flow = allocateFlowRadial(A, heatNode, heatPipe); % 步骤4: 热力计算,供水温度沿管道传播 Ts = propagateSupplyTemp(A, m_flow, Qin, heatPipe); % 步骤5: 根据各节点热负荷计算回水温度 Tr = zeros(nn, 1); for i = 1:nn if heatNode(i).load > 0 && heatNode(i).m_req > 0 Tr(i) = Ts(i) - heatNode(i).load * 1000 / (cp * heatNode(i).m_req); else Tr(i) = Ts(i); end end end实际工程中热网往往是多热源环网,这时候不能简单从单一热源往外推流量,需要求解线性方程组或使用环路平差法。在代码里我做了两种模式的开关:辐射状管网用流量逐级分配,环形管网用回路压降方程联立求解。第一次实现建议先用辐射状,等整个耦合流程跑通之后再扩展环网,这样排查问题更容易定位是热网部分还是耦合迭代部分出错。
3.5 主程序全局耦合迭代逻辑
主程序的耦合迭代是整套代码的灵魂,我把它简化如下:
% 初始化耦合变量 Pchp = 0.3; % CHP电出力 MW Qchp = 0.36; % CHP热出力 MW Peb = 0.4; % 电锅炉耗电功率 MW Qeb = 0.38; % 电锅炉热出力 MW Fgb = 0.1; % 燃气锅炉耗气量 m3/s tolOuter = 1e-4; maxOuter = 30; for k = 1:maxOuter % 保存上一轮耦合功率 Pchp_old = Pchp; Peb_old = Peb; % 更新电网节点注入功率 Ssp(gridNodePqIndex) = -loadGrid + [Pchp; -Peb]; % 求解电网潮流 [V, theta] = solvePowerFlowNR(Ygrid, Ssp, V0, 1e-6, 20); % 更新热网热源功率 QinHeatSource = [Qchp; Qgb; Qeb]; % 求解热网 [Ts, Tr] = solveHeatNetwork(heatNode, heatPipe, QinHeatSource); % 更新气网节点负荷 QgasLoad = [Fchp; Fgb]; % 求解气网 [Pgas] = solveGasNetwork(gasNode, gasPipe, gasSource, QgasLoad); % 根据新的网络状态更新耦合变量 % (“以热定电”模式下CHP电出力受热网回水温度影响) Pchp_new = CHP_ElectricalOutput(Qchp_target, Ts, Tr); Peb_new = ElectricBoilerPower(Qheat_demand, Qchp, Qgb); % 收敛判断 if max(abs(Pchp_new - Pchp), abs(Peb_new - Peb)) < tolOuter break; end Pchp = Pchp_new; Peb = Peb_new; end这个流程的收敛判据只检查CHP和电锅炉的功率变化量,因为它们是连接电网和热网的枢纽设备,只要这两个功率稳定了,整个多能系统就接近稳态了。收敛阈值我取1e-4兆瓦级别,对应热功率0.1千瓦,工程上完全够用。
3.6 天然气网络求解实现注意事项
天然气网的计算对初学者来说容易忽略一个关键点:节点压力基准和气体状态方程的选择。标准Weymouth方程形式为:
Qmn = C·sqrt(Pm² - Pn²)
这里的Qmn是标准状态体积流量,Pm和Pn是绝对压力,且必须以bar为单位才能得到正确的流量量级。如果直接用Pa代入,结果会差好几个数量级。
我在代码里的处理方式是全部用bar作为压力的内部单位,节点压力数据从输入文件读入时自动转换,计算完成后再换算回Pa输出,这样可以避免量纲混乱。另外,天然气网求解初值设置也很重要,我一般把所有节点初压设为气源压力的80%左右,迭代会更容易收敛,这个经验在气网呈辐射状且负荷集中在末端时特别有用。
4. 算例验证与结果分析
4.1 测试系统的具体参数
算例的全部参数我按下面的原则配置。配电网采用IEEE 33节点标准算例,基准电压12.66kV,基准功率10MVA,总负荷约3.7MW,我在此基础上增加了两个耦合设备接入点:节点18接入电锅炉,节点33接入CHP机组。热网是6节点辐射状系统,一个热源节点连接CHP和燃气锅炉,四个负荷节点分别为建筑供暖、工业工艺热、生活热水和农业温室热负荷,总热负荷约1.5MW。气网为6节点系统,两个天然气源节点,一个通向CHP,一个通向燃气锅炉,其余节点作为管道连接点。
CHP参数方面,额定电出力0.3MW,热电比取1.2,发电效率35%,额定热出力0.36MW。电锅炉额定功率0.4MW,效率0.95。燃气锅炉容量0.6MW,效率0.9。热网的供回水设计温度分别为90°C和70°C,水的比热容取4.187kJ/(kg·K),热网节点流量由热负荷除以设计温差得到。
4.2 计算结果展示与收敛行为
程序在18次内层电网迭代和12次外层耦合迭代后收敛,最终结果如下:
- 电网根节点电压0.998pu,节点18电压0.936pu,节点33电压0.952pu,最低电压在节点18为0.931pu,整体电压水平处于可接受范围。
- 节点18的电锅炉从电网吸收有功功率0.38MW,由于电锅炉效率95%,向热网注入0.361MW热功率。
- 节点33的CHP机组以“以热定电”模式运行,电出力0.31MW,热出力0.372MW。
- 热网各负荷节点回水温度最低为64.2°C,最高为68.8°C,均保持在设计允许范围。
- 气网节点最低压力2.22bar,高于供气设备允许的最低进气压,系统整体运行无越限。
迭代曲线显示,前四轮外层迭代中CHP电出力波动幅度较大,从初值0.3MW跳到0.356MW又回落至0.312MW,之后逐渐趋稳。这种初期震荡属于正常现象,在耦合较强的系统中几乎必然出现,因为电网和热网之间存在反馈环路,状态量的调节有一个逐步衰减的过程。
4.3 多能耦合对系统运行状态的实质性影响
这个算例最有价值的结论在于耦合效应不可忽略。我做了对照实验:忽略耦合关系,把电锅炉当成恒定0.4MW纯电网负荷,忽略CHP在电网侧的有功注入,先算电网潮流,再单独算热网。结果显示电网节点18的电压为0.922pu,比耦合计算的结果低0.9个百分点,而CHP节点33的电压则被低估了1.6个百分点。
这个差异在工程上非常明显。电压偏差超过1%就足以影响无功补偿设备的投切策略,更不用说如果电网电压偏低导致电锅炉实际功率下降,热网侧热出力也会跟着减少,这种连锁反应只有耦合计算才能准确揭示。也就是说,在进行综合能源系统的规划设计、设备选型和运行调度时,多能耦合电气热能流计算提供的运行基准数据比分别计算可靠得多。
5. 常见问题与调试经验速查
5.1 外层耦合迭代不收敛怎么办
这是所有跑这种代码的人最先遇到,也最让人头疼的问题。不收敛的原因通常有三个。
第一是耦合变量初值距离真实解太远。解决办法是先用解耦计算给一个粗略初值,也就是先忽略管网影响,直接按能量平衡手动计算一遍CHP和电锅炉的功率,再代入耦合迭代。我实测下来,初值偏差控制在50%以内时收敛基本没问题,偏差超过100%就很容易发散。
第二是热网水力计算结果不更新。常见情况是电锅炉功率上轮算出来是0.38MW,这轮热网算完又要求它输出0.41MW,但代码里热网计算用的还是上一轮的热源数据,导致迭代变量永远追不上目标值。排查方法很简单,在每次外层迭代末尾打印所有耦合设备功率,看看变化量是不是单调递减,如果数值来回跳,多半就是数据更新顺序写错了。
第三是部分状态量越过物理边界导致数值异常,比如节点电压算成负值、热网温度超过200度。遇到这种情况首先要检查迭代步长有没有上限,我建议给状态量的更新加一个阻尼系数,每次只走20%的修正量,虽然收敛慢一点,但稳定性大幅提升。
5.2 雅可比矩阵奇异或条件数过大的排查思路
电网牛拉法中雅可比矩阵奇异通常和节点连接关系、参数设置有关。第一类情况是系统中存在孤立节点,这个节点没有和任何线路相连,导纳矩阵对应行全零,雅可比矩阵必然奇异。解决方法是在结构体数组里检查每条线路的节点编号是否都存在于节点列表中。
第二类情况是PV节点无功越限。配电网中电压支撑能力弱,如果PV节点的无功出力设置超出设备能力,牛拉法迭代中电压修正量会异常,雅可比矩阵逼近奇异。解决方法是加无功越限处理逻辑,一旦PV节点无功越限就把它降级为PQ节点,重新计算。
热网侧的雅可比矩阵问题通常源自关联矩阵秩缺失。当热网存在孤立环时,如果简化成辐射状拓扑强行求解,关联矩阵会出现不满秩。我在代码中增加了拓扑连通性检查,先判断图是否连通,再判断是否存在环,根据结果自动切换求解模式。
5.3 单位换算与数量级引发的错误如何快速定位
综合能源系统最坑人的地方就是单位混用。电网功率用MW和MVar,热网热功率用MW或kW,气网流量用m³/h或kg/s,压力用MPa或bar,温度用°C或K,随便一个换算错误就导致结果差出天价数字。
我的经验是建立一套内部约定的“标准单位制”,代码里一切变量都统一用SI单位加上倍率后缀的命名规则,变量名直接体现单位。比如P_mw表示功率单位是MW,Flow_m3h表示流量单位是m³/h,这样在检查代码时一眼就能发现单位不一致的地方。算例调试时还可以在关键函数入口加断言,判断输入量级是否在合理范围内,比如热功率输入限制在0到10MW之间,超出就报错。这个小检查在项目初期能省下大量排查时间。
5.4 热网供回水温度结果的合理性检查
算完热网后第一件事不是看收敛标志,而是检查温度结果是否符合物理常识。供水温度在热源附近应该最高,沿管道向负荷末端逐渐下降;回水温度在负荷最大节点处最低,然后沿回水管路向热源回升。如果计算结果是供水温度反而低于某个负荷节点的回水温度,说明热网拓扑方向或者流量分配出了问题。
这个问题的常见根源在于管道流量方向设置错误。辐射状热网里供热管道水流方向一定是从热源指向负荷,回水管道方向一定是从负荷指向热源。我在拓扑输入文件中用from和to字段表示管道方向,但一旦有人填反了,关联矩阵的正负号就全错了,计算出的温度传播方向就是颠倒的。
调试时我建议在propagateSupplyTemp函数中打印每个节点的温度值,顺着从热源到末端的方向手动检查一遍,基本能定位到具体是哪根管道方向标反了。
5.5 收敛判据应该怎么选才能既稳又快
收敛判据选得太松,结果失真;选得太严,迭代次数暴增甚至永远无法收敛。电网潮流部分内层收敛阈值我取1e-6功率基准值,也就是10MVA基准下的0.01W,这个精度对应论文出图表完全够用。外层耦合迭代阈值取耦合功率变化量小于1e-4MW,即0.1kW,热功率0.1kW在区域供能系统里大约是整套系统功率的万分之几,作为稳态判据非常合理。
还有一点值得注意:外层迭代的最大次数不宜设置过大,我设为30次,超过即判定不收敛并输出诊断信息。因为交替迭代在正常情况下十几轮就能收敛,如果到30轮还不收敛,继续迭代也只是浪费计算时间,不如停下来查问题。这种“快速失败”的理念在科研代码里非常值得采用,宁可多设置几个断点,也不要让程序默默跑半天最后给你一个错误结果。
6. 项目扩展方向与实际应用思考
这个项目的代码框架有很强的可扩展性。目前实现的是稳态潮流计算,在研究综合能源系统的规划配置、多能互补调度和运行评估时,只需要把外层循环换成时序仿真,对24小时或8760小时的场景逐时段求解,就能得到系统全年的运行轨迹。我在后续工作中就用这套代码做了冬季典型日的逐时仿真,耦合迭代每时段平均收敛时间不到0.5秒,完全可以接受。
另一个扩展方向是加入储能设备和柔性负荷模型。把蓄热电锅炉或蓄冷罐接入热网,在电网侧表现为可平移负荷,在热网侧表现为可转移热源,这种双重身份恰恰需要通过功率更新规则耦合进现有迭代框架。我在扩展时只需要在coupleDeviceModel函数里增加一个储能分支,外层循环在更新功率前先调用储能模型的充放电策略,其余部分基本不用改动。
还有一些细节要提醒:以上所有算法都是在Matlab R2023b环境下调试的,代码没有使用任何工具箱依赖,全部是纯脚本编写。这意味着代码可以在GNU Octave等兼容环境中直接运行,也便于后期翻译成Python或C++做工程部署。如果你要把它嵌入更大的业务系统,建议把网络参数读取部分改为标准化的JSON或CSV格式,方便与外部数据库对接。
对我个人而言,把这块代码从零写通、调平、验证结果,最深的体会是:多能耦合系统本身并不复杂,复杂的是各种模型之间的数据接口。电网和热网各自都有成熟的求解方法,难度全在如何让两个网络的信息顺畅地在耦合设备处交换。只要把耦合关系理清楚、接口格式定义好、收敛判据选得当,整套系统的求解就没有想象中那么吓人。希望这篇博文的经验和代码思路能让你少踩一些我踩过的坑。