压缩空气储能和释能阶段模型➕相关文档文献。

建立了压缩空气储能系统中的压缩机、换热器、储气罐、透平、热水罐等设备的数学模型、 并在 Simulink仿真平台上、 按模块化建模方式完成了系统相关程序编写和仿真模型建立、 包含储能和释能两个阶段的模型

压缩空气储能（CAES）是一个典型的热-流-电多物理场耦合系统。要在 MATLAB/Simulink 中实现完整的储能与释能动态仿真，核心在于对气体状态方程、热力学第一/第二定律以及部件特性曲线的建模。

核心数学模型（理论基础）

在编写代码前，你需要理解以下物理方程，它们是所有模块的基础：

理想气体状态方程（用于储气罐、管道）：
P V = m R T
其中 P 为压力，V 为体积，m 为质量，R 为气体常数（空气取 287 J/kg·K），T 为温度。

质量守恒（用于储气罐）：
frac{dm}{dt} = dot{m{in} - dot{m}{out}

能量守恒（用于储气罐、换热器）：
frac{dU}{dt} = dot{m{in} h{in} - dot{m{out} h{out} + dot{Q}
对于定容储气罐，可推导为：
m c_v frac{dT}{dt} = dot{m}{in} c_p T{in} - dot{m}{out} c_p T{out} + dot{Q}

Simulink 模块化代码实现

下面提供四个核心部件的代码逻辑（使用 MATLAB Function 或 Level-2 MATLAB S-Function）。

压缩机模型 (Compressor)
压缩机消耗电功将空气压缩，产生热量。这里采用基于压比的简化模型。

function [P_out, T_out, W_comp] = fcn_Compressor(P_in, T_in, m_dot, eta_isen)
% 输入:
% P_in: 入口压力 (Pa)
% T_in: 入口温度 (K)
% m_dot: 质量流量 (kg/s)
% eta_isen: 等熵效率
% 输出:
% P_out: 出口压力 (Pa)
% T_out: 出口温度 (K)
% W_comp: 压缩功 (W)

% 假设出口压力由下游决定，或者给定目标压比 beta
beta = 10; % 压比 (示例值)
P_out = P_in * beta;

% 计算等熵出口温度
gamma = 1.4; % 空气比热比
T_isen = T_in * (beta)^((gamma-1)/gamma);

% 计算实际出口温度
T_out = T_in + (T_isen - T_in) / eta_isen;

% 计算压缩功 (轴功)
cp = 1005; % 定压比热容 (J/kg.K)
W_comp = m_dot * cp * (T_out - T_in);

% 注意：实际模型中 P_out 可能需要通过求解器迭代，或者作为输入给定

储气罐模型 (Air Storage Tank)
储气罐是储能的核心，压力和温度随充放气动态变化。这是最复杂的微分方程部分。

function [P, T] = fcn_Tank(m_dot_in, m_dot_out, T_in, T_out, V, R, cv)
% 输入:
% m_dot_in, m_dot_out: 质量流量
% T_in, T_out: 温度
% V: 罐体体积 (m3)
% R, cv: 气体常数和比热
% 输出:
% P, T: 罐内当前压力和温度
% 状态变量: m (质量), U (内能)

persistent m U
if isempty(m)
m = 1; % 初始空气质量 (kg)
U = 200000; % 初始内能 (J)
end

% 质量变化率
dm_dt = m_dot_in - m_dot_out;

% 内能变化率 (假设绝热或考虑热损失 Q_dot)
cp = 1005;
R = 287;
gamma = 1.4;
Q_dot = 0; % 简化：暂不考虑散热

dU_dt = m_dot_in * cp * T_in - m_dot_out * cp * T_out + Q_dot;

% 积分 (欧拉法或使用Simulink Integrator)
m = m + dm_dt * 0.001; % 0.001为仿真步长，实际建议用Integrator模块
U = U + dU_dt * 0.001;

% 计算状态
if m > 0
T = U / (m * cv);
P = m * R * T / V;
else
P = 0;
T = 298;
end

透平模型 (Turbine)
透平在释能阶段工作，将高压空气的热能转化为机械能/电能。

function [P_out, T_out, W_turb] = fcn_Turbine(P_in, T_in, m_dot, eta_isen)
% 输入:
% P_in, T_in: 入口状态
% m_dot: 流量
% eta_isen: 效率
% 输出:
% P_out, T_out: 出口状态
% W_turb: 输出功率

% 假设背压 (出口压力)
P_out = 101325; % 1atm

% 计算等熵膨胀温度
beta = P_in / P_out;
gamma = 1.4;
T_isen = T_in / (beta)^((gamma-1)/gamma);

% 实际出口温度
T_out = T_in - eta_isen * (T_in - T_isen);

% 输出功率
cp = 1005;
W_turb = m_dot * cp * (T_in - T_out);

换热器模型 (Heat Exchanger)
用于冷却压缩热或在释能时预热空气。

function [T_out_cold, T_out_hot] = fcn_HeatExchanger(T_in_cold, T_in_hot, m_dot_cold, m_dot_hot, U, A)
% 简化的 NTU-epsilon 方法
C_min = min(m_dot_cold * 1005, m_dot_hot * 1005);
C_max = max(m_dot_cold * 1005, m_dot_hot * 1005);
C_r = C_min / C_max;

NTU = U * A / C_min; % 传热单元数

% 假设为逆流换热器
epsilon = (1 - exp(-NTU * (1 - C_r))) / (1 - C_r * exp(-NTU * (1 - C_r)));

Q_max = C_min * (T_in_hot - T_in_cold);
Q = epsilon * Q_max;

T_out_cold = T_in_cold + Q / (m_dot_cold * 1005);
T_out_hot = T_in_hot - Q / (m_dot_hot * 1005);

Simulink 模型搭建逻辑（流程图）

要在 Simulink 中实现“储能”和“释能”两个阶段，你需要构建如下逻辑：

主控逻辑：
使用 Switch 模块或 Stateflow 图。
当 Mode == ‘Charge’ 时：
开启压缩机（输入电功率）。
开启换热器（冷却空气）。
空气流向：压缩机 -> 换热器 -> 储气罐。
当 Mode == ‘Discharge’ 时：
开启透平（输出电功率）。
空气流向：储气罐 -> 换热器（加热） -> 透平 -> 排气。

数据接口：
将上述代码放入 MATLAB Function 模块中。
使用 Integrator 模块来积分储气罐的质量和能量方程（替代代码中的 persistent 变量，更符合 Simulink 架构）。

建议的仿真参数初始化

建议初始参数如下：
设备 参数名称 推荐值 备注
储气罐 体积 (V) 1000 m^3 大型CAES

初始压力 50 bar 储能前

最终压力 100 bar 储能后

压缩机 压比 (beta) 5~10 多级压缩

等熵效率 0.85

透平 膨胀比 10

效率 0.88

环境 温度 298 K 25°C

储气罐压力随时间上升 的波形（典型的储能充气阶段），这是压缩空气储能系统中最核心的动态过程。

要复现这个波形，关键在于 储气罐的气体状态方程建模。压力不是线性上升的，但由于储气罐体积大、流量相对稳定，在宏观仿真中常呈现近似线性趋势。

以下是实现该压力变化的核心 MATLAB Function 代码（用于 Simulink）。

储气罐压力动态模型 (MATLAB Function)

这段代码基于 理想气体状态方程 和 质量守恒定律。你需要将该模块放入 Simulink 中，并设置积分器（Integrator）来计算压力变化率。

function dPdt = fcn_AirStorage(Tank_V, m_in, m_out, P_now, T_gas)
% 储气罐压力动态计算
% 输入:
% Tank_V: 储气罐体积 (m^3)
% m_in: 入口质量流量 (kg/s)
% m_out: 出口质量流量 (kg/s) (释能时为正)
% P_now: 当前罐内压力 (Pa)
% T_gas: 罐内气体温度 (K) (假设为绝热或等温过程)
% 输出:
% dPdt: 压力变化率 (Pa/s)

% — 物理常数 —
R_air = 287.1; % 空气气体常数 (J/kg*K)

% — 状态方程推导 —
% 由 PV = mRT 推导
% dm/dt = (P dV/dt + V dP/dt) / (RT) -> 忽略体积变化 dV/dt=0
% dm/dt = (V / (RT)) * dP/dt
% 因此 dP/dt = (RT / V) * dm/dt

% 总质量变化率
dm_dt = m_in - m_out;

% 压力变化率
dPdt = (R_air * T_gas / Tank_V) * dm_dt;

如何在 Simulink 中使用

搭建结构：创建一个子系统，输入端口为 m_in（压缩机来气）、m_out（透平用气）、T_gas（温度，可设为常数如 300K）。
积分器：将上述函数的输出 dPdt 连接到一个 Integrator 模块。
反馈：Integrator 的输出（即当前压力 P）需要反馈回该函数作为 P_now 输入（虽然在等温假设下计算 dPdt 时直接用 dm_dt 即可，但如果涉及变温计算则需要）。
参数设置：
Tank_V：设为 1000 m³（大型储气库）或 10 m³（小型实验台）。
m_in：在储能阶段设为恒定值（如 1 kg/s）。
m_out：在储能阶段设为 0。

换热器（或蓄热罐）温度随时间上升 的过程。温度从初始值（约 100°C）开始，随着时间推移逐渐升高并趋于稳定（约 130°C），这通常发生在压缩空气储能（CAES）系统的 储能阶段，即压缩机排出的高温空气将热量传递给蓄热介质（如水或导热油）。

为了复现这个波形，核心在于 换热器/蓄热罐的热动态模型。这通常是一个一阶惯性环节，描述了热量积累和温度上升的过程。

换热器/蓄热罐温度动态模型 (MATLAB Function)

这段代码基于 热力学能量守恒定律。假设换热器为一个集中参数的热容系统，忽略热损失，或者将其简化为一个一阶惯性环节。

function dTdt = fcn_HeatExchanger(Q_in, T_current, T_amb, C_thermal, R_thermal)
% 换热器/蓄热罐温度动态计算
% 输入:
% Q_in: 输入热量 (W) (来自压缩机排气的热量)
% T_current: 当前换热器温度 (°C 或 K)
% T_amb: 环境温度或冷却介质温度 (°C 或 K) (可选，用于散热)
% C_thermal: 换热器的热容 (J/°C) (质量 * 比热容)
% R_thermal: 热阻 (°C/W) (可选，用于散热到环境)
% 输出:
% dTdt: 温度变化率 (°C/s 或 K/s)

% — 热量平衡方程 —
% 热量积累 = 输入热量 - 散失热量
% C * dT/dt = Q_in - (T - T_amb) / R_thermal

% 如果忽略散热 (R_thermal 设为 inf)，则简化为：
% C * dT/dt = Q_in

% 计算温度变化率
if isnan(R_thermal) || isinf(R_thermal)
% 无散热或散热忽略不计
dTdt = Q_in / C_thermal;
else
% 考虑散热到环境
Q_loss = (T_current - T_amb) / R_thermal;
dTdt = (Q_in - Q_loss) / C_thermal;
end

如何在 Simulink 中使用

搭建结构：创建一个子系统，输入端口为 Q_in（压缩机排气的热量流量）、T_amb（环境温度）、C_thermal（热容）、R_thermal（热阻）。
积分器：将上述函数的输出 dTdt 连接到一个 Integrator 模块，积分器的输出即为换热器的温度 T_current。将 T_current 反馈回函数模块的输入，形成闭环。
参数设置：
C_thermal：根据换热器材料和体积估算。例如，100kg 水的热容约为 100 times 4186 = 4.186 times 10^5 text{ J/°C}。
R_thermal：如果希望温度最终趋于稳定（如你的波形所示），需要设置一个有限的 R_thermal 值，或者在 Q_in 中引入一个随温度变化的项（例如，当温度升高时，换热效率下降，导致 Q_in 减小）。

换热效率下降：随着换热器温度升高，与压缩空气的温差减小，导致传热速率下降。你可以在 Q_in 的计算中引入一个与温差相关的系数。
流量变化：压缩机的排气流量或温度可能不是恒定的，而是随时间变化的。