主题039：污染物排放模型

目录

1. 引言
2. 氮氧化物(NOx)生成机理
3. 一氧化碳(CO)与未燃碳氢化合物
4. 硫氧化物(SOx)生成机理
5. 颗粒物生成模型
6. 污染物控制技术
7. 排放预测模型
8. Python仿真实现
9. 工程应用案例
10. 总结与展望

---

1. 引言

1.1 燃烧污染物的环境影响

燃烧过程是现代社会能源供应的主要方式，但同时也产生了大量的污染物，对环境和人类健康造成了严重影响。主要燃烧污染物包括：

氮氧化物(NOx)：包括NO、NO₂和N₂O，是光化学烟雾、酸雨和臭氧层破坏的主要前驱物。NOx排放会导致：

• 呼吸系统疾病发病率增加
• 光化学烟雾形成，降低能见度
• 酸雨造成土壤和水体酸化
• 温室效应（N₂O的温室效应潜能是CO₂的约300倍）

一氧化碳(CO)：无色无味有毒气体，与血红蛋白的结合能力是氧气的200倍，会导致组织缺氧。主要危害包括：

• 急性中毒甚至死亡
• 慢性暴露导致心血管疾病
• 参与大气化学反应生成臭氧

硫氧化物(SOx)：主要是SO₂和SO₃，是酸雨的主要成分。危害包括：

• 呼吸系统刺激
• 建筑物和文物腐蚀
• 土壤和水体酸化
• 能见度降低

颗粒物(PM)：包括PM₁₀、PM₂.₅和超细颗粒物(UFPs)，可深入人体呼吸系统。危害包括：

• 呼吸系统疾病（哮喘、肺癌）
• 心血管系统疾病
• 能见度降低
• 气候变化（黑碳的辐射强迫效应）

未燃碳氢化合物(UHC)：包括挥发性有机物(VOCs)，是光化学烟雾的前驱物，部分具有致癌性。

1.2 排放法规与标准

随着环保意识的提高，各国制定了严格的排放法规：

中国排放标准：

• 燃煤电厂：NOx ≤ 50 mg/m³，SO₂ ≤ 35 mg/m³，烟尘 ≤ 10 mg/m³（超低排放）
• 燃气轮机：NOx ≤ 50 mg/m³
• 汽车：国六标准，NOx ≤ 0.06 g/km，PM ≤ 0.0045 g/km

欧盟标准：

• 工业排放指令(IED)对各类燃烧设备有详细要求
• 汽车欧7标准将进一步收紧NOx和PM限值

美国标准：

• EPA对各类固定源和移动源有详细规定
• 新源性能标准(NSPS)要求最佳可行控制技术(BACT)

1.3 污染物排放模型的重要性

污染物排放模型在燃烧研究和工程应用中具有重要价值：

1. 设计优化：在设计阶段预测排放水平，优化燃烧器和控制系统
2. 运行优化：指导实际运行参数调整，实现低排放运行
3. 法规遵从：评估是否满足排放法规要求
4. 控制策略评估：评估不同控制技术的效果和经济性
5. 科学研究：深入理解污染物生成机理，开发新型低排放技术

---

2. 氮氧化物(NOx)生成机理

2.1 NOx生成路径概述

燃烧过程中NOx主要通过三种路径生成：

热力型NOx (Thermal NOx)：由空气中的氮气在高温下氧化生成，遵循Zeldovich机理。这是高温燃烧（>1800K）中NOx的主要来源。

快速型NOx (Prompt NOx)：由碳氢化合物燃烧时产生的CH自由基与N₂反应生成，遵循Fenimore机理。在富燃火焰前锋区域显著。

燃料型NOx (Fuel NOx)：由燃料中含有的氮化合物（如煤中的吡啶、吡咯、蛋白质等）转化生成。这是燃煤和生物质燃烧中NOx的主要来源。

2.2 热力型NOx - Zeldovich机理

2.2.1 基元反应

Zeldovich机理包含以下基元反应：

$$\text{O}+{\text{N}}_2\stackrel{k_1}{\to }\text{NO}+\text{N}\Delta H=+75.1\text{ kcal/mol}$$

$$\text{N}+{\text{O}}_2\stackrel{k_2}{\to }\text{NO}+\text{O}\Delta H=-7.7\text{ kcal/mol}$$

$$\text{N}+\text{OH}\stackrel{k_3}{\to }\text{NO}+\text{H}\Delta H=-10.5\text{ kcal/mol}$$

其中，第一个反应是速率控制步骤，具有很高的活化能（约75 kcal/mol），因此热力型NOx生成对温度极其敏感。

2.2.2 反应速率常数

各基元反应的速率常数（Arrhenius形式）：

$$k_1=1.8\times 10^{14}\exp \left(-\frac{38370}{T}\right)[{\text{cm}}^3/(\text{mol}\cdot \text{s})]$$

$$k_2=1.8\times 10^{10}T\exp \left(-\frac{4680}{T}\right)[{\text{cm}}^3/(\text{mol}\cdot \text{s})]$$

$$k_3=7.1\times 10^{13}\exp \left(-\frac{450}{T}\right)[{\text{cm}}^3/(\text{mol}\cdot \text{s})]$$

2.2.3 NO生成速率方程

基于准稳态近似（N原子浓度处于稳态），NO的生成速率为：

$$\frac{d[\text{NO}]}{dt}=2k_1[\text{O}][{\text{N}}_2]\frac{1-\frac{[\text{NO}{]}^2}{K[{\text{O}}_2][{\text{N}}_2]}}{1+\frac{k_{1b}[\text{NO}]}{k_2[{\text{O}}_2]+k_3[\text{OH}]}}$$

其中$K$是平衡常数，$k_{1b}$是反应1的逆反应速率常数。

在初始阶段（NO浓度很低），可以简化为：

$$\frac{d[\text{NO}]}{dt}=2k_1[\text{O}][{\text{N}}_2]$$

2.2.4 扩展Zeldovich机理

更精确的模型还包括N₂O路径和NNH路径：

N₂O路径：
$${\text{N}}_2+\text{O}+\text{M}\to {\text{N}}_2\text{O}+\text{M}$$
$${\text{N}}_2\text{O}+\text{O}\to 2\text{NO}$$

此路径在中等温度（1500-1800K）和贫燃条件下较重要。

NNH路径：
$${\text{N}}_2+\text{H}\to \text{NNH}$$
$$\text{NNH}+\text{O}\to \text{NO}+\text{NH}$$

此路径在富燃条件下有一定贡献。

2.3 快速型NOx - Fenimore机理

2.3.1 反应机理

快速型NOx由Fenimore于1971年发现，主要反应路径为：

$$\text{CH}+{\text{N}}_2\to \text{HCN}+\text{N}$$
$${\text{CH}}_2+{\text{N}}_2\to \text{HCN}+\text{NH}$$

生成的HCN和NH随后通过一系列反应转化为NO：

$$\text{HCN}+\text{O}\to \text{NCO}+\text{H}$$
$$\text{NCO}+\text{H}\to \text{NH}+\text{CO}$$
$$\text{NH}+\text{H}\to \text{N}+{\text{H}}_2$$
$$\text{N}+\text{OH}\to \text{NO}+\text{H}$$

2.3.2 影响因素

快速型NOx生成受以下因素影响：

1. 当量比：在富燃条件（φ > 1.1）下显著增加
2. 燃料类型：含碳量高的燃料生成更多快速型NOx
3. 火焰结构：层流火焰前锋区域贡献最大
4. 温度：虽然活化能较低，但仍随温度升高而增加

快速型NOx通常占总NOx的10-30%，但在某些条件下（如富燃天然气燃烧）可达50%以上。

2.4 燃料型NOx

2.4.1 燃料氮的赋存形态

燃料中的氮主要以以下形态存在：

煤中的氮：

• 吡咯型氮（N-5）：五元环结构，约占40-60%
• 吡啶型氮（N-6）：六元环结构，约占20-40%
• 季氮（N-Q）：约占5-20%
• 胺型氮：少量

生物质中的氮：

• 蛋白质氮
• 氨基酸
• 叶绿素

重油中的氮：

• 碱性氮（吡啶类）
• 非碱性氮（吡咯类、咔唑类）

2.4.2 转化机理

燃料氮向NOx的转化包括以下步骤：

1. 挥发分氮释放：燃料热解析出含氮挥发分（HCN、NH₃等）
2. 挥发分氮氧化：挥发分中的氮化合物氧化生成NO
3. 焦炭氮氧化：焦炭中的氮在燃烧过程中转化为NO
4. NO还原：在还原性气氛中，NO可被还原为N₂

2.4.3 转化模型

燃料型NOx生成可以用以下简化模型描述：

$$\text{燃料N}\stackrel{\beta }{\to }\text{挥发分N}+\text{焦炭N}$$

$$\text{挥发分N}\stackrel{\gamma }{\to }\text{NO}\text{(挥发分N转化率)}$$

$$\text{焦炭N}\stackrel{\delta }{\to }\text{NO}\text{(焦炭N转化率)}$$

总燃料型NOx生成：

$$[\text{NO}{]}_{\text{fuel}}=\beta \cdot \gamma \cdot N_{\text{vol}}+(1-\beta )\cdot \delta \cdot N_{\text{char}}$$

其中：

• β：挥发分氮比例（通常0.6-0.8）
• γ：挥发分氮转化率（通常0.2-0.5）
• δ：焦炭氮转化率（通常0.1-0.3）

转化率受温度、氧浓度、停留时间等因素影响。

2.5 NOx生成经验模型

2.5.1 热力型NOx经验公式

基于大量实验数据，可以建立热力型NOx生成的经验关联式：

指数型关联式：

$$[{\text{NO}}_x]=A\cdot \exp \left(-\frac{E_a}{RT}\right)\cdot [{\text{O}}_2{]}^n\cdot {\tau }^m$$

其中：

• A：前置因子
• Ea：表观活化能（约60-70 kcal/mol）
• n：氧浓度指数（通常0.5-1.0）
• m：停留时间指数（通常0.5-1.0）
• τ：停留时间

De Soete关联式：

$$\frac{d[\text{NO}]}{dt}=1.35\times 10^{14}\exp \left(-\frac{67500}{RT}\right)[{\text{O}}_2{]}^{0.5}[{\text{N}}_2]$$

2.5.2 燃料型NOx经验公式

燃料型NOx可以用以下关联式估算：

$${\eta }_{\text{NOx}}=\frac{[{\text{NO}}_x]}{N_{\text{fuel}}}=f(T,\varphi ,\tau ,N_{\text{fuel}})$$

其中η_NOx是燃料氮向NOx的转化率，通常在10-40%范围内。

2.6 NOx还原机理

2.6.1 热还原（Reburning）

在还原性气氛中，NO可以被还原：

$$\text{NO}+{\text{CH}}_i\to \text{HCN}+\text{产物}$$
$$\text{NO}+{\text{NH}}_i\to {\text{N}}_2+\text{产物}$$
$$\text{NO}+\text{CO}\to \frac{1}{2}{\text{N}}_2+{\text{CO}}_2$$

再燃技术就是利用这一原理，在主燃烧区后喷入再燃燃料，创造还原性气氛还原已生成的NOx。

2.6.2 选择性非催化还原(SNCR)

使用氨或尿素作为还原剂，在无催化剂条件下选择性还原NOx：

$$4\text{NO}+4{\text{NH}}_3+{\text{O}}_2\to 4{\text{N}}_2+6{\text{H}}_2\text{O}$$
$$2{\text{NO}}_2+4{\text{NH}}_3+{\text{O}}_2\to 3{\text{N}}_2+6{\text{H}}_2\text{O}$$

最佳反应温度窗口：850-1100°C。温度过低反应速率慢，温度过高氨会被氧化生成NOx。

2.6.3 选择性催化还原(SCR)

在催化剂作用下，氨选择性还原NOx：

$$4\text{NO}+4{\text{NH}}_3+{\text{O}}_2\stackrel{\text{催化剂}}{\to }4{\text{N}}_2+6{\text{H}}_2\text{O}$$

常用催化剂：V₂O₅-WO₃/TiO₂、沸石等。

反应温度窗口：300-400°C，脱硝效率可达90%以上。

---

3. 一氧化碳(CO)与未燃碳氢化合物

3.1 CO生成机理

3.1.1 化学动力学路径

CO是碳氢燃料氧化的中间产物，主要生成路径包括：

甲醛路径：
$${\text{CH}}_2\text{O}+\text{OH}\to \text{HCO}+{\text{H}}_2\text{O}$$
$$\text{HCO}+\text{M}\to \text{H}+\text{CO}+\text{M}$$

乙醛路径：
$${\text{CH}}_3\text{CHO}+\text{OH}\to {\text{CH}}_3\text{CO}+{\text{H}}_2\text{O}$$
$${\text{CH}}_3\text{CO}+\text{M}\to {\text{CH}}_3+\text{CO}+\text{M}$$

直接氧化：
$$\text{C}+\frac{1}{2}{\text{O}}_2\to \text{CO}$$

3.1.2 CO氧化为CO₂

CO最终氧化为CO₂的反应：

$$\text{CO}+\text{OH}\to {\text{CO}}_2+\text{H}$$

这是CO氧化的主要路径，反应速率：

$$\frac{d[\text{CO}]}{dt}=-k[\text{CO}][\text{OH}]$$

其中：

$$k=2.74\times 10^9{T}^{1.3}\exp \left(-\frac{1160}{T}\right)[{\text{cm}}^3/(\text{mol}\cdot \text{s})]$$

3.1.3 CO生成的影响因素

1. 温度：低温下CO氧化速率慢，CO排放高
2. 混合程度：混合不良导致局部富燃，CO生成增加
3. 停留时间：停留时间不足，CO来不及完全氧化
4. 当量比：富燃条件下CO排放显著增加

3.2 未燃碳氢化合物(UHC)

3.2.1 UHC来源

未燃碳氢化合物的主要来源：

1. 淬熄效应：火焰接触冷壁面或冷流体而熄灭，燃料未完全燃烧
2. 混合不良：燃料与氧化剂混合不充分，形成富燃区域
3. 缝隙效应：火焰无法进入狭窄缝隙，燃料逃逸
4. 不完全燃烧：氧气不足或温度过低导致燃烧不完全

3.2.2 UHC排放特征

UHC排放特征：

• 冷启动时UHC排放最高
• 怠速和低负荷时UHC排放较高
• 过量空气系数偏离化学计量比时UHC增加
• 燃烧室设计不良导致UHC升高

3.3 CO与UHC的关联

CO和UHC排放通常同时出现，但机理不同：

• CO：燃烧中间产物，氧化不完全导致
• UHC：燃料未参与燃烧或燃烧不完全

在贫燃条件下，CO通常先于UHC被氧化；在富燃条件下，两者同时增加。

---

4. 硫氧化物(SOx)生成机理

4.1 燃料硫的形态

燃料中的硫以多种形态存在：

煤中的硫：

• 有机硫（与碳结合的硫，如噻吩、硫醇）
• 黄铁矿硫（FeS₂）
• 硫酸盐硫（CaSO₄等，通常不燃烧）

石油中的硫：

• 硫醇（R-SH）
• 硫醚（R-S-R）
• 噻吩及其衍生物
• 二硫化物（R-S-S-R）

天然气中的硫：

• 主要是H₂S和有机硫化合物

4.2 SOx生成机理

4.2.1 SO₂生成

燃料硫燃烧主要生成SO₂：

$$\text{S}+{\text{O}}_2\to {\text{SO}}_2$$
$${\text{H}}_2\text{S}+\frac{3}{2}{\text{O}}_2\to {\text{SO}}_2+{\text{H}}_2\text{O}$$
$$\text{COS}+\frac{3}{2}{\text{O}}_2\to {\text{SO}}_2+{\text{CO}}_2$$

SO₂生成量与燃料硫含量成正比：

$$[{\text{SO}}_2]=2\times \frac{S_{\text{fuel}}}{M_S}\times \frac{M_{{\text{SO}}_2}}{V_{\text{flue}}}$$

其中：

• S_fuel：燃料硫含量（质量分数）
• M_S：硫的摩尔质量（32 g/mol）
• M_SO₂：SO₂的摩尔质量（64 g/mol）
• V_flue：单位质量燃料产生的烟气体积

4.2.2 SO₃生成

少量SO₂进一步氧化生成SO₃：

$${\text{SO}}_2+\frac{1}{2}{\text{O}}_2\to {\text{SO}}_3$$

SO₃生成机理：

1. 气相氧化：高温下SO₂与O原子或OH自由基反应
2. 催化氧化：金属氧化物（如V₂O₅）催化SO₂氧化

SO₃通常占SOx的1-5%，但危害更大：

• 与水反应生成硫酸雾（H₂SO₄），腐蚀设备
• 形成蓝烟，造成可见污染
• 是酸雨的主要成分

4.2.3 SO₃生成动力学

SO₃生成速率可以用以下方程描述：

$$\frac{d[{\text{SO}}_3]}{dt}=k[{\text{SO}}_2][{\text{O}}_2{]}^{0.5}[\text{M}]$$

其中k是速率常数，M是第三体。

SO₃/SO₂平衡比与温度关系：

$$K_p=\frac{[{\text{SO}}_3]}{[{\text{SO}}_2][{\text{O}}_2{]}^{0.5}}$$

温度越低，平衡越偏向SO₃。但在高温燃烧区，平衡偏向SO₂。

4.3 SOx控制技术

4.3.1 燃烧前脱硫

• 煤炭洗选：物理方法去除部分黄铁矿硫
• 煤气化：将煤转化为合成气，同时脱除H₂S
• 煤液化：将煤转化为液体燃料，脱硫更彻底

4.3.2 燃烧中脱硫

炉内喷钙：向炉膛内喷入石灰石（CaCO₃）或生石灰（CaO）

$${\text{CaCO}}_3\to \text{CaO}+{\text{CO}}_2$$
$$\text{CaO}+{\text{SO}}_2+\frac{1}{2}{\text{O}}_2\to {\text{CaSO}}_4$$

脱硫效率：30-50%（循环流化床可达90%以上）

4.3.3 烟气脱硫(FGD)

湿法脱硫（最常用）：

• 石灰石-石膏法：
  $${\text{CaCO}}_3+{\text{SO}}_2+\frac{1}{2}{\text{O}}_2+2{\text{H}}_2\text{O}\to {\text{CaSO}}_4\cdot 2{\text{H}}_2\text{O}+{\text{CO}}_2$$
• 氨法：
  $$2{\text{NH}}_3+{\text{SO}}_2+{\text{H}}_2\text{O}\to ({\text{NH}}_4{)}_2{\text{SO}}_3$$

脱硫效率：95%以上

干法/半干法脱硫：

• 喷雾干燥法
• 循环流化床脱硫
• 电子束脱硫

---

5. 颗粒物生成模型

5.1 颗粒物分类

燃烧产生的颗粒物按粒径分类：

总悬浮颗粒物(TSP)：粒径 < 100 μm
可吸入颗粒物(PM₁₀)：粒径 < 10 μm，可进入呼吸道
细颗粒物(PM₂.₅)：粒径 < 2.5 μm，可进入肺泡
超细颗粒物(UFP)：粒径 < 0.1 μm，可进入血液循环

5.2 颗粒物来源

5.2.1 灰分颗粒

煤和生物质燃烧产生的无机灰分颗粒：

• 飞灰：随烟气排出的细小灰分颗粒
• 底灰：炉膛底部排出的大颗粒灰分

灰分颗粒的形成过程：

1. 矿物颗粒在高温下熔融或部分熔融
2. 挥发分释放导致颗粒破碎
3. 熔融颗粒碰撞聚并
4. 烟气冷却过程中凝固形成玻璃质颗粒

5.2.2 碳烟(Soot)

碳烟是碳氢燃料不完全燃烧产生的碳质颗粒，主要成分为碳（>90%）和少量氢。

碳烟生成机理：

1. 成核：多环芳烃(PAH)通过碰撞和化学反应形成初始颗粒
2. 表面生长：碳烟颗粒表面通过乙炔(C₂H₂)等碳源添加碳原子
3. 聚并：小颗粒碰撞聚并形成大颗粒
4. 氧化：碳烟颗粒被O₂、OH等氧化剂氧化消耗

碳烟生成模型：

经验模型（如Khan-Greeves模型）：

$$\frac{d{Y}_s}{dt}=A\cdot p^{0.5}\cdot \exp \left(-\frac{E_a}{RT}\right)\cdot [\text{fuel}{]}^n-B\cdot Y_s\cdot \frac{p_{{\text{O}}_2}}{p}$$

其中：

• Y_s：碳烟质量分数
• A, B：模型常数
• p：压力
• p_O₂：氧分压

详细模型（如PAH生长模型）：

• 跟踪PAH物种浓度
• 考虑PAH成核、生长、氧化全过程
• 计算量大，但更准确

5.2.3 硫酸盐和重金属

• 硫酸盐：SO₃与碱金属反应生成硫酸盐颗粒
• 重金属：燃料中的重金属（Hg、Pb、As等）挥发后冷凝形成颗粒

5.3 颗粒物形成动力学

5.3.1 成核

成核速率可以用经典成核理论描述：

$$J=N_1\sqrt{\frac{2\sigma }{\pi m_1}}\exp \left(-\frac{16\pi {\sigma }^3{v}_1^2}{3(k_{B}T{)}^3{\ln }^{2}S}\right)$$

其中：

• J：成核速率
• N₁：单体浓度
• σ：表面张力
• m₁：单体质量
• v₁：单体体积
• S：过饱和度

5.3.2 表面生长

表面生长速率：

$$\frac{dm}{dt}=\pi d_p^2\cdot {\rho }_s\cdot k_s\cdot [{\text{C}}_2{\text{H}}_2]$$

其中：

• d_p：颗粒直径
• ρ_s：碳烟密度
• k_s：表面生长速率常数

5.3.3 聚并

聚并导致颗粒数减少、粒径增大：

$$\frac{dn}{dt}=-\frac{1}{2}\beta n^2$$

其中：

• n：颗粒数浓度
• β：聚并系数

5.3.4 氧化

碳烟氧化速率：

$$\frac{dm}{dt}=-\pi d_p^2\cdot {\rho }_s\cdot k_{ox}\cdot [{\text{O}}_2]\text{ 或 }[\text{OH}]$$

常用氧化速率表达式（Nagle-Strickland-Constable模型）：

$$k_{ox}=\left(\frac{k_A{p}_{{\text{O}}_2}}{1+k_Z{p}_{{\text{O}}_2}}\right)x+k_B{p}_{{\text{O}}_2}(1-x)$$

其中x是活性表面位点比例。

5.4 粒径分布模型

5.4.1 对数正态分布

颗粒物粒径通常符合对数正态分布：

$$n(d_p)=\frac{N}{\sqrt{2\pi }{d}_p\ln {\sigma }_g}\exp \left(-\frac{(\ln d_p-\ln d_{pg}{)}^2}{2{\ln }^2{\sigma }_g}\right)$$

其中：

• N：总颗粒数
• d_pg：几何平均直径
• σ_g：几何标准偏差

5.4.2 矩方法

矩方法通过求解颗粒分布的各阶矩来简化计算：

$$M_k={\int }_0^{\infty }{d}_p^{k}n(d_p)d{d}_p$$

常用的是M₀（颗粒数）、M₁（总直径）、M₂（总表面积）、M₃（总体积）。

5.5 颗粒物控制技术

5.5.1 机械除尘

• 重力除尘器：利用重力沉降，效率低，用于预处理
• 惯性除尘器：利用惯性分离，效率中等
• 旋风除尘器：利用离心力分离，效率80-90%

5.5.2 电除尘(ESP)

利用高压电场使颗粒带电，然后在集尘极被捕集：

• 除尘效率：99%以上
• 可捕集PM₂.₅和亚微米颗粒
• 能耗相对较低

5.5.3 袋式除尘

利用滤袋过滤烟气中的颗粒：

• 除尘效率：99.9%以上
• 对细颗粒捕集效果好
• 运行阻力较大

5.5.4 湿法除尘

利用液体（通常是水）洗涤烟气：

• 可同时去除颗粒物和SO₂
• 适用于高温高湿烟气
• 产生废水需要处理

---

6. 污染物控制技术

6.1 低NOx燃烧技术

6.1.1 空气分级燃烧

将燃烧所需空气分两级或多级送入：

• 主燃区：富燃燃烧，抑制NOx生成
• 燃尽区：补充空气，完成燃烧

NOx降低幅度：30-50%

6.1.2 燃料分级燃烧（再燃）

在主燃区后喷入再燃燃料（通常是天然气或煤粉）：

• 创造还原性气氛
• 还原已生成的NOx
• 再燃燃料比例：10-25%

NOx降低幅度：50-70%

6.1.3 烟气再循环(FGR)

将部分烟气返回燃烧区：

• 降低火焰温度，抑制热力型NOx
• 降低氧浓度，减缓NOx生成

NOx降低幅度：20-50%

6.1.4 低NOx燃烧器

通过特殊设计实现低NOx燃烧：

• 旋流燃烧器：强化混合，缩短高温停留时间
• 分级燃烧器：实现空气或燃料分级
• 富燃-淬熄-贫燃( RQL)燃烧器：航空发动机常用

6.2 后处理技术

6.2.1 选择性催化还原(SCR)

在催化剂作用下，氨选择性还原NOx：

• 反应温度：300-400°C
• 脱硝效率：90%以上
• 氨逃逸：❤️ ppm

催化剂类型：

• 钒基催化剂（V₂O₅-WO₃/TiO₂）
• 沸石催化剂（ZSM-5、β-沸石）
• 贵金属催化剂（Pt、Pd）

6.2.2 选择性非催化还原(SNCR)

在无催化剂条件下，氨或尿素还原NOx：

• 反应温度：850-1100°C
• 脱硝效率：30-50%
• 投资成本低于SCR

6.2.3 联合脱硝技术

SNCR-SCR联合：

• SNCR作为预处理，降低SCR负荷
• 减少催化剂用量
• 提高整体脱硝效率

低NOx燃烧+SCR：

• 低NOx燃烧降低初始NOx浓度
• SCR深度脱硝
• 实现超低排放（NOx < 50 mg/m³）

6.3 多污染物协同控制

6.3.1 活性炭/焦吸附

活性炭或活性焦可同时吸附多种污染物：

• 脱除SO₂、NOx、Hg、二噁英等
• 再生后可重复使用
• 适用于小型锅炉和垃圾焚烧

6.3.2 臭氧氧化

利用臭氧将NO氧化为易溶的高价氮氧化物：

$$\text{NO}+{\text{O}}_3\to {\text{NO}}_2+{\text{O}}_2$$
$${\text{NO}}_2+{\text{O}}_3\to {\text{N}}_2{\text{O}}_5+{\text{O}}_2$$

然后在湿法脱硫系统中脱除。

6.3.3 电子束法

利用电子束辐射烟气，产生自由基氧化污染物：

• 同时脱除SO₂和NOx
• 副产物可作肥料
• 能耗较高，应用受限

---

7. 排放预测模型

7.1 零维模型

7.1.1 化学平衡模型

假设燃烧达到化学平衡，计算平衡产物组成：

$$\min G=\sum _i{n}_i{\mu }_i$$

约束条件：

• 元素守恒
• 质量守恒
• 能量守恒

适用于估算理论最低排放水平。

7.1.2 反应器模型

充分搅拌反应器(PSR)：

$$\frac{d{n}_i}{dt}={\dot{n}}_{i,\text{in}}-{\dot{n}}_{i,\text{out}}+V\cdot {\dot{\omega }}_i=0$$

适用于模拟燃烧室中的充分混合区域。

柱塞流反应器(PFR)：

$$\frac{d{n}_i}{dx}=A\cdot {\dot{\omega }}_i$$

适用于模拟一维流动和反应。

7.2 一维模型

7.2.1 火焰模型

求解一维层流火焰方程：

连续性方程：
$$\dot{m}=\rho uA=\text{const}$$

能量方程：
$$\dot{m}\frac{dT}{dx}-\frac{d}{dx}\left(\lambda \frac{dT}{dx}\right)+\sum _i\rho Y_i{V}_i{c}_{p,i}\frac{dT}{dx}+\sum _i{h}_i{\dot{\omega }}_i{W}_i=0$$

组分方程：
$$\dot{m}\frac{d{Y}_i}{dx}+\frac{d}{dx}(\rho Y_i{V}_i)-{\dot{\omega }}_i{W}_i=0$$

可预测火焰中NOx、CO等污染物的生成分布。

7.3 CFD耦合模型

7.3.1 湍流燃烧模型

涡耗散模型(EDM)：

$${\dot{\omega }}_{\text{fuel}}=A\rho \frac{\epsilon }{k}\min \left(Y_{\text{fuel}},\frac{Y_{\text{ox}}}{s}\right)$$

假设反应速率由混合控制，适用于快速化学反应。

火焰面模型：

在混合分数空间求解标量输运：
$$\rho \frac{\partial \varphi }{\partial t}+\rho u\frac{\partial \varphi }{\partial x}=\frac{\partial }{\partial x}\left(\rho D\frac{\partial \varphi }{\partial x}\right)+{\dot{\omega }}_{\varphi }$$

适用于非预混火焰。

7.3.2 污染物后处理模型

NOx后处理模型：

在流场计算完成后，基于温度和组分场计算NOx生成：

$$\frac{d[\text{NO}]}{dt}=f(T,[{\text{O}}_2],[{\text{N}}_2],\tau )$$

采用"冻结流"假设，即NO生成不影响主燃烧过程。

碳烟模型：

求解碳烟输运方程：

$$\frac{\partial \rho Y_s}{\partial t}+\nabla \cdot (\rho u{Y}_s)=\nabla \cdot (\rho D\nabla Y_s)+S_s$$

其中S_s是碳烟源项（成核+生长-氧化）。

7.4 机器学习模型

7.4.1 代理模型

使用神经网络等机器学习方法建立排放预测代理模型：

输入参数：

• 燃料性质（热值、组分）
• 运行参数（负荷、空燃比）
• 设计参数（燃烧器类型、尺寸）

输出：

• NOx、CO、SOx、颗粒物排放浓度

优势：

• 计算速度快（毫秒级）
• 可用于实时优化
• 可嵌入控制系统

7.4.2 模型训练

训练数据来源于：

• 实验测量
• 详细CFD计算
• 文献数据

常用算法：

• 人工神经网络(ANN)
• 支持向量机(SVM)
• 随机森林(RF)
• 高斯过程回归(GPR)

---

8. Python仿真实现

8.1 热力型NOx计算

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint

# 设置中文字体
plt.rcParams['font.size'] = 10

# 物理常数
R = 8.314  # J/(mol·K)

# Zeldovich机理速率常数 (cm^3/(mol·s))
def k1(T):
    """O + N2 -> NO + N"""
    return 1.8e14 * np.exp(-38370 / T)

def k2(T):
    """N + O2 -> NO + O"""
    return 1.8e10 * T * np.exp(-4680 / T)

def k3(T):
    """N + OH -> NO + H"""
    return 7.1e13 * np.exp(-450 / T)

def k1b(T):
    """NO + N -> O + N2 (逆反应)"""
    return 3.8e12 * np.exp(-425 / T)

# 平衡常数
def K_NO(T):
    """NO生成反应平衡常数"""
    return 4.0 * np.exp(-21650 / T)

# 热力型NOx生成速率
def thermal_NO_rate(T, O2, N2, NO, O, OH):
    """
    计算热力型NOx生成速率
    单位: mol/(cm^3·s)
    """
    # 使用准稳态近似
    # 简化模型：d[NO]/dt = 2*k1*[O]*[N2]
    
    # O原子平衡浓度（来自O2解离）
    O_eq = 3.0e16 * O2**0.5 * np.exp(-27100 / T)  # 近似表达式
    
    # NO生成速率
    if NO < 1e-10:  # 初始阶段
        dNO_dt = 2 * k1(T) * O_eq * N2 * 1e-6  # 转换为mol/(cm^3·s)
    else:
        # 考虑逆反应
        K = K_NO(T)
        dNO_dt = 2 * k1(T) * O_eq * N2 * (1 - NO**2 / (K * O2 * N2)) * 1e-6
    
    return max(0, dNO_dt)

# 计算NOppm随时间变化
def calculate_NO_history(T, O2_frac, N2_frac, P=1.0, t_max=1.0):
    """
    计算NO浓度随时间变化
    T: 温度 (K)
    O2_frac: O2摩尔分数
    N2_frac: N2摩尔分数
    P: 压力 (atm)
    t_max: 最大时间 (s)
    """
    # 总浓度 (mol/cm^3)
    C_total = P * 1.013e5 / (R * T) * 1e-6  # mol/cm^3
    
    O2 = O2_frac * C_total
    N2 = N2_frac * C_total
    
    # 时间数组
    t = np.linspace(0, t_max, 1000)
    
    # 求解ODE
    def dNO_dt(NO, t):
        O_eq = 3.0e16 * O2**0.5 * np.exp(-27100 / T) * 1e-6  # mol/cm^3
        rate = 2 * k1(T) * O_eq * N2
        
        # 考虑平衡限制
        K = K_NO(T)
        if NO > 0 and O2 > 0 and N2 > 0:
            rate *= max(0, 1 - NO**2 / (K * O2 * N2))
        
        return rate
    
    NO = odeint(dNO_dt, 0, t)
    
    # 转换为ppm
    NO_ppm = NO / C_total * 1e6
    
    return t, NO_ppm.flatten()

# 实验1: 温度对热力型NOx的影响
print("=" * 60)
print("实验1: 温度对热力型NOx生成的影响")
print("=" * 60)

temperatures = np.linspace(1500, 2500, 11)
NO_final = []

fig, axes = plt.subplots(2, 2, figsize=(12, 10))

# 计算不同温度下的NO生成
for T in temperatures[::2]:  # 每隔一个温度计算
    t, NO_ppm = calculate_NO_history(T, 0.05, 0.75, t_max=2.0)
    axes[0, 0].semilogy(t, NO_ppm, label=f'T={T:.0f}K')
    NO_final.append(NO_ppm[-1])

axes[0, 0].set_xlabel('Time (s)')
axes[0, 0].set_ylabel('NO Concentration (ppm)')
axes[0, 0].set_title('Thermal NO Formation at Different Temperatures')
axes[0, 0].legend()
axes[0, 0].grid(True, alpha=0.3)

# 最终NO浓度vs温度
axes[0, 1].semilogy(temperatures, 
                    [calculate_NO_history(T, 0.05, 0.75, t_max=2.0)[1][-1] 
                     for T in temperatures], 
                    'ro-', linewidth=2, markersize=8)
axes[0, 1].axvline(x=1800, color='g', linestyle='--', label='Typical threshold (1800K)')
axes[0, 1].set_xlabel('Temperature (K)')
axes[0, 1].set_ylabel('Final NO Concentration (ppm)')
axes[0, 1].set_title('Temperature Effect on Thermal NOx')
axes[0, 1].legend()
axes[0, 1].grid(True, alpha=0.3)

print("\n温度对热力型NOx的影响:")
for T in [1600, 1800, 2000, 2200, 2400]:
    _, NO_ppm = calculate_NO_history(T, 0.05, 0.75, t_max=2.0)
    print(f"  T = {T}K: NO = {NO_ppm[-1]:.2f} ppm")

plt.tight_layout()
plt.savefig('exp1_thermal_NO_temperature.png', dpi=150, bbox_inches='tight')
print("\n结果已保存至: exp1_thermal_NO_temperature.png")
plt.close()

8.2 燃料型NOx计算

# 实验2: 燃料型NOx生成模型
print("\n" + "=" * 60)
print("实验2: 燃料型NOx生成模型")
print("=" * 60)

def fuel_NOx_model(N_fuel, T, phi, tau, beta=0.7, gamma=0.35, delta=0.2):
    """
    燃料型NOx生成模型
    
    参数:
    N_fuel: 燃料氮含量 (wt%)
    T: 温度 (K)
    phi: 当量比
    tau: 停留时间 (s)
    beta: 挥发分氮比例
    gamma: 挥发分氮转化率
    delta: 焦炭氮转化率
    
    返回:
    NOx_ppm: NOx浓度 (ppm, 干基6%O2)
    """
    # 温度修正因子
    f_T = np.exp(-(T - 1500)**2 / (2 * 200**2)) + 0.3
    
    # 当量比修正因子 (富燃时转化率降低)
    if phi > 1.0:
        f_phi = np.exp(-0.5 * (phi - 1.0))
    else:
        f_phi = 1.0 - 0.2 * (1.0 - phi)
    
    # 停留时间修正因子
    f_tau = 1.0 - np.exp(-tau / 0.5)
    
    # 总转化率
    eta_vol = gamma * f_T * f_phi * f_tau
    eta_char = delta * f_T * f_phi * f_tau
    
    # 燃料氮向NOx的总转化率
    eta_total = beta * eta_vol + (1 - beta) * eta_char
    
    # 假设燃料完全燃烧，计算烟气量
    # 简化：1kg燃料产生约10Nm³烟气
    V_flue = 10.0  # Nm³/kg燃料
    
    # 燃料氮质量 (kg/kg燃料)
    m_N = N_fuel / 100.0
    
    # NOx质量 (以NO₂计)
    m_NOx = m_N * eta_total * (46/14)  # NO₂分子量/氮原子量
    
    # 转换为mg/Nm³
    NOx_mg = m_NOx * 1e6 / V_flue  # mg/Nm³
    
    # 转换为ppm (干基6%O2)
    # 假设烟气中NOx浓度与O2浓度关系
    NOx_ppm = NOx_mg / 2.05  # 近似转换
    
    return NOx_ppm, eta_total

# 计算不同条件下的燃料型NOx
N_fuel_values = [0.5, 1.0, 1.5, 2.0, 2.5]  # 燃料氮含量 (%)
T_values = np.linspace(1200, 1800, 20)  # 温度范围
phi_values = np.linspace(0.8, 1.4, 15)  # 当量比范围

fig, axes = plt.subplots(2, 2, figsize=(12, 10))

# 燃料氮含量影响
NOx_vs_N = []
for N in N_fuel_values:
    NOx, _ = fuel_NOx_model(N, 1500, 1.0, 1.0)
    NOx_vs_N.append(NOx)
    print(f"  燃料氮含量 {N}%: NOx = {NOx:.1f} ppm")

axes[0, 0].plot(N_fuel_values, NOx_vs_N, 'bo-', linewidth=2, markersize=10)
axes[0, 0].set_xlabel('Fuel Nitrogen Content (wt%)')
axes[0, 0].set_ylabel('NOx Emission (ppm)')
axes[0, 0].set_title('Effect of Fuel Nitrogen Content on NOx')
axes[0, 0].grid(True, alpha=0.3)

# 温度影响
NOx_vs_T = []
eta_vs_T = []
for T in T_values:
    NOx, eta = fuel_NOx_model(1.5, T, 1.0, 1.0)
    NOx_vs_T.append(NOx)
    eta_vs_T.append(eta * 100)

axes[0, 1].plot(T_values, NOx_vs_T, 'r-', linewidth=2, label='NOx')
ax_twin = axes[0, 1].twinx()
ax_twin.plot(T_values, eta_vs_T, 'g--', linewidth=2, label='Conversion Rate')
axes[0, 1].set_xlabel('Temperature (K)')
axes[0, 1].set_ylabel('NOx Emission (ppm)', color='r')
ax_twin.set_ylabel('Conversion Rate (%)', color='g')
axes[0, 1].set_title('Effect of Temperature on Fuel NOx')
axes[0, 1].grid(True, alpha=0.3)

# 当量比影响
NOx_vs_phi = []
for phi in phi_values:
    NOx, _ = fuel_NOx_model(1.5, 1500, phi, 1.0)
    NOx_vs_phi.append(NOx)

axes[1, 0].plot(phi_values, NOx_vs_phi, 'mo-', linewidth=2, markersize=6)
axes[1, 0].axvline(x=1.0, color='k', linestyle='--', label='Stoichiometric')
axes[1, 0].set_xlabel('Equivalence Ratio')
axes[1, 0].set_ylabel('NOx Emission (ppm)')
axes[1, 0].set_title('Effect of Equivalence Ratio on Fuel NOx')
axes[1, 0].legend()
axes[1, 0].grid(True, alpha=0.3)

# 停留时间影响
tau_values = np.linspace(0.1, 3.0, 20)
NOx_vs_tau = []
for tau in tau_values:
    NOx, _ = fuel_NOx_model(1.5, 1500, 1.0, tau)
    NOx_vs_tau.append(NOx)

axes[1, 1].plot(tau_values, NOx_vs_tau, 'c-', linewidth=2)
axes[1, 1].set_xlabel('Residence Time (s)')
axes[1, 1].set_ylabel('NOx Emission (ppm)')
axes[1, 1].set_title('Effect of Residence Time on Fuel NOx')
axes[1, 1].grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig('exp2_fuel_NOx_model.png', dpi=150, bbox_inches='tight')
print("\n结果已保存至: exp2_fuel_NOx_model.png")
plt.close()

8.3 碳烟生成模型

# 实验3: 碳烟生成模型
print("\n" + "=" * 60)
print("实验3: 碳烟生成模型")
print("=" * 60)

def soot_model(T, phi, C2H2, O2, OH, dt=0.001, t_max=0.1):
    """
    简化碳烟生成模型
    
    参数:
    T: 温度 (K)
    phi: 当量比
    C2H2: 乙炔浓度 (mol/cm^3)
    O2: O2浓度 (mol/cm^3)
    OH: OH浓度 (mol/cm^3)
    
    返回:
    t: 时间数组
    Y_soot: 碳烟质量分数历史
    n_soot: 颗粒数浓度历史
    d_p: 颗粒直径历史
    """
    # 模型参数
    A_nuc = 1.0e5  # 成核前置因子
    Ea_nuc = 15000  # 成核活化能 (cal/mol)
    A_sg = 1.0e4   # 表面生长前置因子
    Ea_sg = 12000  # 表面生长活化能
    A_coag = 1.0e-10  # 聚并系数
    A_ox = 1.0e3   # 氧化前置因子
    Ea_ox = 18000  # 氧化活化能
    
    # 碳烟密度 (g/cm^3)
    rho_s = 1.8
    
    # 碳烟分子量 (g/mol)
    W_s = 12.0
    
    # 初始条件
    Y_soot = 0.0  # 碳烟质量分数
    n_soot = 1.0e6  # 颗粒数浓度 (#/cm^3)
    
    # 时间积分
    t = np.arange(0, t_max, dt)
    Y_history = []
    n_history = []
    d_history = []
    
    for ti in t:
        # 成核速率 (#/cm^3/s)
        k_nuc = A_nuc * np.exp(-Ea_nuc / (1.987 * T))
        J_nuc = k_nuc * C2H2**2
        
        # 表面生长速率 (g/cm^3/s)
        k_sg = A_sg * np.exp(-Ea_sg / (1.987 * T))
        
        # 聚并速率
        beta_coag = A_coag * np.sqrt(T)
        
        # 氧化速率 (g/cm^3/s)
        k_ox = A_ox * np.exp(-Ea_ox / (1.987 * T))
        
        # 计算颗粒直径
        if n_soot > 1.0e3 and Y_soot > 1.0e-15:
            d_p = (6 * Y_soot / (np.pi * rho_s * n_soot))**(1/3)
        else:
            d_p = 1.0e-7  # 初始直径 1 nm
        
        # 碳烟质量变化
        if Y_soot > 0 and n_soot > 1.0e3:
            surface_area = np.pi * d_p**2 * n_soot
            dY_dt = (k_sg * surface_area * C2H2 - 
                    k_ox * surface_area * (O2 + OH*10))
        else:
            dY_dt = 0
        
        # 颗粒数变化
        dn_dt = J_nuc - 0.5 * beta_coag * n_soot**2
        
        # 更新
        Y_soot = max(0, Y_soot + dY_dt * dt)
        n_soot = max(1.0e3, n_soot + dn_dt * dt)
        
        Y_history.append(Y_soot)
        n_history.append(n_soot)
        d_history.append(d_p * 1e7)  # 转换为nm
    
    return t, np.array(Y_history), np.array(n_history), np.array(d_history)

# 计算不同当量比下的碳烟生成
phi_values = [1.5, 1.8, 2.0, 2.2, 2.5]
T = 1800  # K

fig, axes = plt.subplots(2, 2, figsize=(12, 10))

colors = plt.cm.viridis(np.linspace(0, 1, len(phi_values)))

for i, phi in enumerate(phi_values):
    # 假设乙炔浓度随当量比增加
    C2H2 = 1.0e-8 * phi**2
    O2 = 1.0e-7 * max(0, 2.0 - phi) / phi
    OH = 1.0e-9 * np.exp(-(phi - 1.0)**2)
    
    t, Y_soot, n_soot, d_p = soot_model(T, phi, C2H2, O2, OH, t_max=0.05)
    
    axes[0, 0].semilogy(t*1000, Y_soot*1e6, color=colors[i], 
                        label=f'φ={phi}', linewidth=2)
    axes[0, 1].semilogy(t*1000, n_soot, color=colors[i], 
                        label=f'φ={phi}', linewidth=2)
    axes[1, 0].plot(t*1000, d_p, color=colors[i], 
                    label=f'φ={phi}', linewidth=2)

axes[0, 0].set_xlabel('Time (ms)')
axes[0, 0].set_ylabel('Soot Mass Fraction (ppm)')
axes[0, 0].set_title('Soot Mass Fraction Evolution')
axes[0, 0].legend()
axes[0, 0].grid(True, alpha=0.3)

axes[0, 1].set_xlabel('Time (ms)')
axes[0, 1].set_ylabel('Soot Number Density (#/cm³)')
axes[0, 1].set_title('Soot Number Density Evolution')
axes[0, 1].legend()
axes[0, 1].grid(True, alpha=0.3)

axes[1, 0].set_xlabel('Time (ms)')
axes[1, 0].set_ylabel('Particle Diameter (nm)')
axes[1, 0].set_title('Soot Particle Growth')
axes[1, 0].legend()
axes[1, 0].grid(True, alpha=0.3)

# 温度对碳烟生成的影响
T_values = np.linspace(1400, 2200, 20)
Y_soot_final = []

for T in T_values:
    phi = 2.0
    C2H2 = 1.0e-8 * phi**2
    O2 = 1.0e-7 * max(0, 2.0 - phi) / phi
    OH = 1.0e-9 * np.exp(-(phi - 1.0)**2)
    
    t, Y_soot, n_soot, d_p = soot_model(T, phi, C2H2, O2, OH, t_max=0.05)
    Y_soot_final.append(Y_soot[-1] * 1e6)

axes[1, 1].plot(T_values, Y_soot_final, 'r-', linewidth=2)
axes[1, 1].set_xlabel('Temperature (K)')
axes[1, 1].set_ylabel('Final Soot Mass Fraction (ppm)')
axes[1, 1].set_title('Temperature Effect on Soot Formation')
axes[1, 1].grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig('exp3_soot_formation.png', dpi=150, bbox_inches='tight')
print("\n结果已保存至: exp3_soot_formation.png")
plt.close()

print("\n碳烟生成模型结果:")
print(f"  富燃条件(φ=2.0, T=1800K)下最终碳烟浓度: {Y_soot_final[8]:.2f} ppm")

8.4 综合排放预测

# 实验4: 综合排放预测模型
print("\n" + "=" * 60)
print("实验4: 燃烧污染物排放综合预测")
print("=" * 60)

def comprehensive_emission_model(fuel_type, load, excess_air, T_combustor, 
                                  N_content=0, S_content=0):
    """
    综合燃烧污染物排放预测模型
    
    参数:
    fuel_type: 燃料类型 ('natural_gas', 'coal', 'diesel')
    load: 负荷率 (0-1)
    excess_air: 过量空气系数
    T_combustor: 燃烧室温度 (K)
    N_content: 燃料氮含量 (wt%)
    S_content: 燃料硫含量 (wt%)
    
    返回:
    emissions: 污染物排放字典 (ppm或mg/m³)
    """
    
    # 当量比
    phi = 1.0 / excess_air
    
    # O2浓度 (假设)
    O2_conc = 0.21 * (excess_air - 1) / excess_air
    
    # N2浓度
    N2_conc = 0.79
    
    emissions = {}
    
    # ========== NOx预测 ==========
    # 热力型NOx
    if T_combustor > 1700:
        thermal_NOx = 50 * np.exp((T_combustor - 1800) / 100) * load
    else:
        thermal_NOx = 5 * load
    
    # 快速型NOx (富燃时增加)
    if phi > 1.0:
        prompt_NOx = 20 * (phi - 1.0) * load
    else:
        prompt_NOx = 5 * load
    
    # 燃料型NOx
    if N_content > 0:
        fuel_NOx = N_content * 50 * load  # 简化模型
    else:
        fuel_NOx = 0
    
    emissions['NOx'] = thermal_NOx + prompt_NOx + fuel_NOx
    
    # ========== CO预测 ==========
    # 低温或富燃时CO增加
    if T_combustor < 1500:
        CO_factor = 1500 / T_combustor
    else:
        CO_factor = 1.0
    
    if phi > 1.1:
        CO_rich = 500 * (phi - 1.0)
    else:
        CO_rich = 10
    
    emissions['CO'] = CO_rich * CO_factor * (1.5 - load)  # 低负荷时CO高
    
    # ========== SOx预测 ==========
    if S_content > 0:
        emissions['SO2'] = S_content * 20 * load  # mg/m³
    else:
        emissions['SO2'] = 0
    
    # ========== 颗粒物预测 ==========
    if fuel_type == 'coal':
        base_PM = 50 * load
    elif fuel_type == 'diesel':
        base_PM = 30 * load if phi > 1.2 else 10 * load
    else:
        base_PM = 5 * load
    
    # 富燃增加碳烟
    if phi > 1.2:
        soot_PM = 100 * (phi - 1.2)**2 * load
    else:
        soot_PM = 0
    
    emissions['PM'] = base_PM + soot_PM  # mg/m³
    
    # ========== UHC预测 ==========
    if phi > 1.1 or T_combustor < 1400:
        emissions['UHC'] = 50 * (phi - 0.9) * (1500 / max(T_combustor, 1200))
    else:
        emissions['UHC'] = 5
    
    return emissions

# 模拟不同燃烧条件
conditions = [
    {'name': '天然气锅炉-满负荷', 'fuel': 'natural_gas', 'load': 1.0, 
     'air': 1.2, 'T': 1900, 'N': 0, 'S': 0},
    {'name': '天然气锅炉-半负荷', 'fuel': 'natural_gas', 'load': 0.5, 
     'air': 1.3, 'T': 1700, 'N': 0, 'S': 0},
    {'name': '燃煤锅炉-满负荷', 'fuel': 'coal', 'load': 1.0, 
     'air': 1.25, 'T': 1650, 'N': 1.2, 'S': 2.5},
    {'name': '柴油机-高负荷', 'fuel': 'diesel', 'load': 1.0, 
     'air': 1.5, 'T': 2000, 'N': 0.1, 'S': 0.3},
    {'name': '柴油机-低负荷', 'fuel': 'diesel', 'load': 0.3, 
     'air': 2.0, 'T': 1500, 'N': 0.1, 'S': 0.3},
]

print("\n不同燃烧条件下的污染物排放预测:")
print("-" * 80)
print(f"{'工况':<20} {'NOx(ppm)':<12} {'CO(ppm)':<12} {'SO2(mg/m³)':<15} {'PM(mg/m³)':<12} {'UHC(ppm)':<10}")
print("-" * 80)

results = []
for cond in conditions:
    emis = comprehensive_emission_model(
        cond['fuel'], cond['load'], cond['air'], 
        cond['T'], cond['N'], cond['S']
    )
    results.append({
        'name': cond['name'],
        'NOx': emis['NOx'],
        'CO': emis['CO'],
        'SO2': emis['SO2'],
        'PM': emis['PM'],
        'UHC': emis['UHC']
    })
    print(f"{cond['name']:<20} {emis['NOx']:<12.1f} {emis['CO']:<12.1f} "
          f"{emis['SO2']:<15.1f} {emis['PM']:<12.1f} {emis['UHC']:<10.1f}")

# 可视化
fig, axes = plt.subplots(2, 3, figsize=(15, 10))

names = [r['name'] for r in results]
x_pos = np.arange(len(names))

# NOx
axes[0, 0].bar(x_pos, [r['NOx'] for r in results], color='red', alpha=0.7)
axes[0, 0].set_xticks(x_pos)
axes[0, 0].set_xticklabels(names, rotation=45, ha='right')
axes[0, 0].set_ylabel('NOx (ppm)')
axes[0, 0].set_title('NOx Emissions')
axes[0, 0].grid(True, alpha=0.3, axis='y')

# CO
axes[0, 1].bar(x_pos, [r['CO'] for r in results], color='orange', alpha=0.7)
axes[0, 1].set_xticks(x_pos)
axes[0, 1].set_xticklabels(names, rotation=45, ha='right')
axes[0, 1].set_ylabel('CO (ppm)')
axes[0, 1].set_title('CO Emissions')
axes[0, 1].grid(True, alpha=0.3, axis='y')

# SO2
axes[0, 2].bar(x_pos, [r['SO2'] for r in results], color='blue', alpha=0.7)
axes[0, 2].set_xticks(x_pos)
axes[0, 2].set_xticklabels(names, rotation=45, ha='right')
axes[0, 2].set_ylabel('SO2 (mg/m³)')
axes[0, 2].set_title('SO2 Emissions')
axes[0, 2].grid(True, alpha=0.3, axis='y')

# PM
axes[1, 0].bar(x_pos, [r['PM'] for r in results], color='gray', alpha=0.7)
axes[1, 0].set_xticks(x_pos)
axes[1, 0].set_xticklabels(names, rotation=45, ha='right')
axes[1, 0].set_ylabel('PM (mg/m³)')
axes[1, 0].set_title('Particulate Matter Emissions')
axes[1, 0].grid(True, alpha=0.3, axis='y')

# UHC
axes[1, 1].bar(x_pos, [r['UHC'] for r in results], color='green', alpha=0.7)
axes[1, 1].set_xticks(x_pos)
axes[1, 1].set_xticklabels(names, rotation=45, ha='right')
axes[1, 1].set_ylabel('UHC (ppm)')
axes[1, 1].set_title('Unburned HC Emissions')
axes[1, 1].grid(True, alpha=0.3, axis='y')

# 雷达图 - 综合排放比较
ax_radar = axes[1, 2]
ax_radar.remove()
ax_radar = fig.add_subplot(2, 3, 6, projection='polar')

categories = ['NOx', 'CO', 'SO2', 'PM', 'UHC']
N = len(categories)

# 归一化数据
angles = [n / float(N) * 2 * np.pi for n in range(N)]
angles += angles[:1]

colors_radar = ['red', 'blue', 'green', 'orange', 'purple']
for i, cond in enumerate(conditions[:3]):  # 只画前3个
    emis = comprehensive_emission_model(
        cond['fuel'], cond['load'], cond['air'], 
        cond['T'], cond['N'], cond['S']
    )
    values = [
        emis['NOx'] / 200,  # 归一化
        emis['CO'] / 100,
        emis['SO2'] / 100,
        emis['PM'] / 100,
        emis['UHC'] / 50
    ]
    values += values[:1]
    
    ax_radar.plot(angles, values, 'o-', linewidth=2, 
                  label=cond['name'], color=colors_radar[i])
    ax_radar.fill(angles, values, alpha=0.25, color=colors_radar[i])

ax_radar.set_xticks(angles[:-1])
ax_radar.set_xticklabels(categories)
ax_radar.set_title('Emission Comparison (Normalized)')
ax_radar.legend(loc='upper right', bbox_to_anchor=(1.3, 1.0))

plt.tight_layout()
plt.savefig('exp4_comprehensive_emissions.png', dpi=150, bbox_inches='tight')
print("\n结果已保存至: exp4_comprehensive_emissions.png")
plt.close()

print("\n综合排放预测结果分析:")
print("  1. 天然气锅炉NOx主要来自热力型，温度控制是关键")
print("  2. 燃煤锅炉燃料型NOx占主导，需考虑低氮燃烧+SCR")
print("  3. 柴油机PM排放显著，需要颗粒物捕集器")
print("  4. 低负荷运行CO和UHC排放增加")

8.5 控制技术效果评估

# 实验5: 污染物控制技术效果评估
print("\n" + "=" * 60)
print("实验5: 污染物控制技术效果评估")
print("=" * 60)

def evaluate_control_technologies(baseline_emissions, technologies):
    """
    评估不同控制技术的效果
    
    参数:
    baseline_emissions: 基准排放字典
    technologies: 控制技术列表，每项包含名称和脱除效率
    
    返回:
    控制后的排放字典
    """
    results = {'Baseline': baseline_emissions.copy()}
    
    current_emissions = baseline_emissions.copy()
    
    for tech in technologies:
        name = tech['name']
        efficiency = tech['efficiency']  # 各污染物脱除效率字典
        
        controlled = {}
        for pollutant in current_emissions:
            if pollutant in efficiency:
                removal = efficiency[pollutant]
                controlled[pollutant] = current_emissions[pollutant] * (1 - removal)
            else:
                controlled[pollutant] = current_emissions[pollutant]
        
        results[name] = controlled.copy()
        current_emissions = controlled.copy()
    
    return results

# 燃煤锅炉基准排放
baseline = {
    'NOx': 400,  # ppm
    'SO2': 1500,  # mg/m³
    'PM': 80,   # mg/m³
    'CO': 50    # ppm
}

# 控制技术序列
technologies = [
    {
        'name': 'Low-NOx Burner',
        'efficiency': {'NOx': 0.30}
    },
    {
        'name': '+ SNCR',
        'efficiency': {'NOx': 0.40}
    },
    {
        'name': '+ SCR',
        'efficiency': {'NOx': 0.85}
    },
    {
        'name': '+ FGD',
        'efficiency': {'SO2': 0.95}
    },
    {
        'name': '+ ESP',
        'efficiency': {'PM': 0.995}
    }
]

control_results = evaluate_control_technologies(baseline, technologies)

# 可视化
fig, axes = plt.subplots(2, 2, figsize=(12, 10))

stages = list(control_results.keys())
x_pos = np.arange(len(stages))

# NOx控制效果
NOx_values = [control_results[s]['NOx'] for s in stages]
axes[0, 0].bar(x_pos, NOx_values, color='red', alpha=0.7)
axes[0, 0].axhline(y=50, color='g', linestyle='--', label='Ultra-low emission limit')
axes[0, 0].set_xticks(x_pos)
axes[0, 0].set_xticklabels(stages, rotation=45, ha='right')
axes[0, 0].set_ylabel('NOx (ppm)')
axes[0, 0].set_title('NOx Control Technology Roadmap')
axes[0, 0].legend()
axes[0, 0].grid(True, alpha=0.3, axis='y')

# SO2控制效果
SO2_values = [control_results[s]['SO2'] for s in stages]
axes[0, 1].bar(x_pos, SO2_values, color='blue', alpha=0.7)
axes[0, 1].axhline(y=35, color='g', linestyle='--', label='Ultra-low emission limit')
axes[0, 1].set_xticks(x_pos)
axes[0, 1].set_xticklabels(stages, rotation=45, ha='right')
axes[0, 1].set_ylabel('SO2 (mg/m³)')
axes[0, 1].set_title('SO2 Control Technology Roadmap')
axes[0, 1].legend()
axes[0, 1].grid(True, alpha=0.3, axis='y')

# PM控制效果
PM_values = [control_results[s]['PM'] for s in stages]
axes[1, 0].bar(x_pos, PM_values, color='gray', alpha=0.7)
axes[1, 0].axhline(y=10, color='g', linestyle='--', label='Ultra-low emission limit')
axes[1, 0].set_xticks(x_pos)
axes[1, 0].set_xticklabels(stages, rotation=45, ha='right')
axes[1, 0].set_ylabel('PM (mg/m³)')
axes[1, 0].set_title('PM Control Technology Roadmap')
axes[1, 0].legend()
axes[1, 0].grid(True, alpha=0.3, axis='y')

# 综合脱除效率
pollutants = ['NOx', 'SO2', 'PM']
removal_efficiency = []
for p in pollutants:
    initial = control_results['Baseline'][p]
    final = control_results[stages[-1]][p]
    removal = (initial - final) / initial * 100
    removal_efficiency.append(removal)

axes[1, 1].bar(pollutants, removal_efficiency, color=['red', 'blue', 'gray'], alpha=0.7)
axes[1, 1].set_ylabel('Removal Efficiency (%)')
axes[1, 1].set_title('Overall Pollutant Removal Efficiency')
axes[1, 1].grid(True, alpha=0.3, axis='y')

for i, v in enumerate(removal_efficiency):
    axes[1, 1].text(i, v + 1, f'{v:.1f}%', ha='center', va='bottom', fontweight='bold')

plt.tight_layout()
plt.savefig('exp5_control_technologies.png', dpi=150, bbox_inches='tight')
print("\n结果已保存至: exp5_control_technologies.png")
plt.close()

print("\n控制技术效果评估:")
print("-" * 60)
for stage in stages:
    print(f"\n{stage}:")
    for pollutant, value in control_results[stage].items():
        print(f"  {pollutant}: {value:.1f}")

8.6 排放法规符合性分析

# 实验6: 排放法规符合性分析
print("\n" + "=" * 60)
print("实验6: 排放法规符合性分析")
print("=" * 60)

# 中国超低排放标准 (燃煤电厂)
ultra_low_standards = {
    'NOx': 50,    # mg/m³
    'SO2': 35,    # mg/m³
    'PM': 10,     # mg/m³
    'Hg': 0.03    # mg/m³
}

# 不同工况下的排放
operating_conditions = [
    {'name': '100% Load', 'NOx': 45, 'SO2': 30, 'PM': 8, 'Hg': 0.025},
    {'name': '75% Load', 'NOx': 48, 'SO2': 32, 'PM': 9, 'Hg': 0.028},
    {'name': '50% Load', 'NOx': 52, 'SO2': 35, 'PM': 10, 'Hg': 0.030},
    {'name': '30% Load', 'NOx': 58, 'SO2': 38, 'PM': 12, 'Hg': 0.035},
]

# 单位转换 (ppm -> mg/m³ for NOx)
# NOx: 1 ppm ≈ 2.05 mg/m³ (标准状态)
for cond in operating_conditions:
    cond['NOx_mg'] = cond['NOx'] * 2.05

fig, axes = plt.subplots(2, 2, figsize=(12, 10))

conditions_names = [c['name'] for c in operating_conditions]
x_pos = np.arange(len(conditions_names))

# NOx符合性
NOx_values = [c['NOx_mg'] for c in operating_conditions]
axes[0, 0].bar(x_pos, NOx_values, color='red', alpha=0.7)
axes[0, 0].axhline(y=ultra_low_standards['NOx'], color='g', 
                   linestyle='--', linewidth=2, label='Standard (50 mg/m³)')
axes[0, 0].set_xticks(x_pos)
axes[0, 0].set_xticklabels(conditions_names)
axes[0, 0].set_ylabel('NOx (mg/m³)')
axes[0, 0].set_title('NOx Emission Compliance')
axes[0, 0].legend()
axes[0, 0].grid(True, alpha=0.3, axis='y')

# 标记超标
for i, v in enumerate(NOx_values):
    if v > ultra_low_standards['NOx']:
        axes[0, 0].text(i, v + 1, 'Non-compliant', ha='center', 
                        color='red', fontsize=8, rotation=45)

# SO2符合性
SO2_values = [c['SO2'] for c in operating_conditions]
axes[0, 1].bar(x_pos, SO2_values, color='blue', alpha=0.7)
axes[0, 1].axhline(y=ultra_low_standards['SO2'], color='g', 
                   linestyle='--', linewidth=2, label='Standard (35 mg/m³)')
axes[0, 1].set_xticks(x_pos)
axes[0, 1].set_xticklabels(conditions_names)
axes[0, 1].set_ylabel('SO2 (mg/m³)')
axes[0, 1].set_title('SO2 Emission Compliance')
axes[0, 1].legend()
axes[0, 1].grid(True, alpha=0.3, axis='y')

# PM符合性
PM_values = [c['PM'] for c in operating_conditions]
axes[1, 0].bar(x_pos, PM_values, color='gray', alpha=0.7)
axes[1, 0].axhline(y=ultra_low_standards['PM'], color='g', 
                   linestyle='--', linewidth=2, label='Standard (10 mg/m³)')
axes[1, 0].set_xticks(x_pos)
axes[1, 0].set_xticklabels(conditions_names)
axes[1, 0].set_ylabel('PM (mg/m³)')
axes[1, 0].set_title('PM Emission Compliance')
axes[1, 0].legend()
axes[1, 0].grid(True, alpha=0.3, axis='y')

# 综合符合性指数
compliance_index = []
for cond in operating_conditions:
    idx = 0
    idx += min(1.0, ultra_low_standards['NOx'] / cond['NOx_mg'])
    idx += min(1.0, ultra_low_standards['SO2'] / cond['SO2'])
    idx += min(1.0, ultra_low_standards['PM'] / cond['PM'])
    idx += min(1.0, ultra_low_standards['Hg'] / cond['Hg'])
    compliance_index.append(idx / 4 * 100)  # 百分比

colors_comp = ['green' if c >= 100 else 'orange' if c >= 90 else 'red' 
               for c in compliance_index]
axes[1, 1].bar(x_pos, compliance_index, color=colors_comp, alpha=0.7)
axes[1, 1].axhline(y=100, color='k', linestyle='-', linewidth=2, label='Full Compliance')
axes[1, 1].axhline(y=90, color='orange', linestyle='--', label='Warning Level')
axes[1, 1].set_xticks(x_pos)
axes[1, 1].set_xticklabels(conditions_names)
axes[1, 1].set_ylabel('Compliance Index (%)')
axes[1, 1].set_title('Overall Emission Compliance Index')
axes[1, 1].set_ylim([0, 120])
axes[1, 1].legend()
axes[1, 1].grid(True, alpha=0.3, axis='y')

plt.tight_layout()
plt.savefig('exp6_emission_compliance.png', dpi=150, bbox_inches='tight')
print("\n结果已保存至: exp6_emission_compliance.png")
plt.close()

print("\n排放法规符合性分析:")
print("-" * 60)
print(f"{'工况':<15} {'NOx':<12} {'SO2':<12} {'PM':<12} {'符合性':<10}")
print("-" * 60)
for i, cond in enumerate(operating_conditions):
    status = '✓' if compliance_index[i] >= 100 else '✗'
    print(f"{cond['name']:<15} {cond['NOx_mg']:<12.1f} {cond['SO2']:<12.1f} "
          f"{cond['PM']:<12.1f} {status:<10}")

---

9. 工程应用案例

9.1 燃煤电厂超低排放改造

9.1.1 项目背景

某600MW燃煤电厂进行超低排放改造，目标达到：

• NOx ≤ 50 mg/m³
• SO₂ ≤ 35 mg/m³
• 烟尘 ≤ 10 mg/m³

9.1.2 改造方案

低NOx燃烧系统：

• 采用低NOx燃烧器，降低热力型NOx生成
• 空气分级燃烧，主燃区富燃运行
• 烟气再循环，降低火焰温度

SNCR+SCR联合脱硝：

• SNCR作为预处理，在850-1000°C区域喷入氨水
• SCR深度脱硝，使用V₂O₅-WO₃/TiO₂催化剂
• 脱硝效率达到90%以上

湿法脱硫：

• 石灰石-石膏湿法脱硫
• 脱硫效率95%以上
• 增设GGH消除白烟

除尘系统：

• 低低温电除尘器+湿式电除尘器
• 除尘效率99.9%以上

9.1.3 改造效果

改造后排放水平：

• NOx: 35-45 mg/m³
• SO₂: 20-30 mg/m³
• 烟尘: 3-8 mg/m³

达到超低排放要求，年减少污染物排放：

• NOx: 约800吨
• SO₂: 约3000吨
• 烟尘: 约400吨

9.2 燃气轮机低NOx燃烧器设计

9.2.1 技术挑战

燃气轮机燃烧室温度高达2000K以上，热力型NOx生成显著。主要挑战：

• 高温燃烧导致NOx排放高
• 燃烧稳定性与低排放的平衡
• 变工况运行的适应性

9.2.2 设计方案

干式低NOx(DLN)燃烧器：

• 燃料分级：主燃料+预混燃料
• 贫燃预混燃烧，降低火焰温度
• 旋流稳定火焰

燃烧模式：

• 启动：扩散燃烧模式，保证点火和稳定性
• 低负荷：部分预混模式
• 高负荷：完全预混模式，最低NOx排放

9.2.3 性能指标

• NOx排放: < 25 ppm (@15% O₂)
• CO排放: < 10 ppm
• 燃烧效率: > 99.9%
• 燃烧脉动: < 1% (RMS)

9.3 柴油机排放后处理系统

9.3.1 系统组成

氧化催化器(DOC)：

• 氧化CO和HC
• 将NO氧化为NO₂（促进DPF再生）
• 贵金属催化剂(Pt/Pd)

颗粒物捕集器(DPF)：

• 壁流式蜂窝陶瓷
• 捕集效率>95%
• 定期再生（被动+主动）

选择性催化还原(SCR)：

• 尿素喷射系统
• 钒基或沸石催化剂
• NOx转化效率>90%

氨逃逸催化器(ASC)：

• 氧化逃逸的NH₃
• 防止二次污染

9.3.2 国六排放标准达标

后处理系统使柴油机满足国六b标准：

• NOx: 0.06 g/kWh (相比国五降低80%)
• PM: 0.01 g/kWh (相比国五降低50%)
• PN: 6×10¹¹ #/kWh

---