梯度校正法燃气热值测量装置抽气系统建模

刘 伟，王玉刚，张洪军，赵晓东，潘 江

(中国计量大学计量测试工程学院，浙江 杭州 310018)

0 引 言

在天然气能量计量方法中，以燃烧法测量可燃气体热值是能源领域的基础问题，也是目前研究热点。国外很多学者和计量机构已经在这方面有较为深入的研究并取得长足进步。20世纪30年代美国国家标准局的Rossini建立了电能标定的等温式热值测量计，以电校正的形式将热值溯源到质量、温度、电流和时间等国际单位上[1-3]。20世纪80年代，中国计量科学研究院在Rossini方法的基础上用氧弹法测量纯甲烷的热值，于2012年开始研制基准气体热量计[4-5]。此后，中国石油西南油气田分公司天然气研究院、中国计量大学等研究机构或高校，在天然气热值测量技术、配气比例、热值测量主体、热量计非稳态温度场优化分析、抽气系统设计等方面取得较大进展[6-10]。中国计量大学研制的气体热值测量系统是在Rossini法的基础上改进而来，目前该装置已初步成型。但该系统在抽气控制方面的研究还尚待完善，如只能实现恒速抽气，没有建立完善的数学模型，无法应对诸如外界气压变化、管路内部气源不稳、抽气过快或过慢等不利因素引起的燃烧不稳定情况。本文的研究内容就是建立抽气系统的准确数学模型，并对模型进行验证，为后期设计数字控制器从而实现抽气系统变速抽气控制做准备。

1 抽气系统工作原理

抽气系统结构如图1所示，系统由干燥管、缓冲瓶、尾气收集装置、压差传感器、数据采集单元、单片机电路以及控制程序组成。系统工作时，操作人员通过控制程序设置好抽气装置运行参数，发布抽气命令，单片机接收到参数和命令后，控制尾气收集装置运行的相关电路导通，尾气收集装置按照设定的运行速度将燃烧室排气口出来的尾气经过干燥后抽入尾气收集装置的气缸内；压差传感器检测到燃烧室进、排气口两端压差的变化并输出电压，使用Agilent 34972A数据采集单元采集该电压，通过串口将电压发送到控制程序中再换算为气体压差值，经过控制器的计算后给出下一时刻的运行速度信号，再通过串口自动向单片机发布新的运行速度信号；如此循环往复，在控制程序的控制下实现运行速度的自调整，在一次燃烧实验完成后，排出气缸内的尾气。

图1 抽气系统结构示意图

2 实验信号设计

2.1 M序列生成

为了预先对抽气系统的动态特性有一定了解，如确定抽气系统是否存在滞后现象，确定系统调节过程时间等，对抽气系统进行阶跃响应测试[11]。阶跃响应测试时，抽气系统工作在开环控制状态。设定尾气收集装置活塞运行速度为20 mm/min，待抽气系统运行处于稳定抽气状态时，保持当前运行速度，在此运行速度基础上再突加一个幅值为β的正向速度信号，作为阶跃信号输入抽气系统。为确定多大幅值的β值能获得较为理想的阶跃响应信号，在20 mm/min运行速度基础上分别施加原有运行速度的5%、10%、15%、20%、25%、50%作为阶跃信号输入系统，抽气系统在不同幅值阶跃信号测试下的响应曲线如图2所示。实验结果表明，当阶跃信号幅值β选取为原来速度的25%，即β=5 mm/min时，系统阶跃响应效果较为理想。

图2 不同幅值阶跃信号测试下系统的响应曲线

在得到不同幅值阶跃信号测试结果后，将β的幅值设定为5 mm/min再次进行测试，记录数据如图3所示。多次重复实验结果表明，抽气系统的调节时间Ts≈21 s，纯滞后时间 τ0≈2 s。

图3 5 mm/min阶跃信号测试下系统的响应曲线

冲激信号作为最理想的系统测试信号，拥有广泛的频谱，但让计算机产生一个理想的冲激信号较为困难。研究发现，白噪声的自相关函数与冲激信号很接近。白噪声取代冲激信号对一个系统进行测试是最为理想的，但在实际工程应用中，执行器并不可能完全按照白噪声的变化规律来动作，而M序列（即最长线性移位寄存器序列）能较好地模拟执行器动态运行中的输入信号。通过对M序列自相关分析发现，M序列自相关函数也接近于冲激信号，其自相关函数图形分析如图4所示。

图4 M序列自相关函数图形

使用M序列产生一系列‘0’和‘1’的码值，并给码值‘0’和码值‘1’赋以相应的活塞运行速度值作为抽气系统激励信号输入系统，得到一系列M序列的响应值。利用移位寄存器产生‘0’和‘1’的码值，需要确定M序列的时钟周期T0、序列长度Np、序列幅值A0和试验长度N。已知抽气系统的调节时间Ts，M序列的时钟周期可根据式（1）来确定：

序列长度根据式（2）确定：

代入数据计算后取M序列时钟周期T0=0.5 s，序列长度Np=63。依据Np=26–1=63，确定M序列由6级移位寄存器发生，6级移位寄存器连接方式用式（3）表示，寄存器连接示意图如图5所示。

图5 6级移位寄存器连接示意图

根据抽气系统阶跃响应情况可知，系统初始运行时会有较长的调节时间，在获取动态数据时，需连续发送两个M序列实验信号。第一个M序列用来消除调节过程，得到的系统响应数据不用于参数辨识，系统的第二个M序列响应数据用于参数辨识。依据6级移位寄存器结构示意图，利用LabVIEW编写M序列码生成程序[12]，初始序列设为101010，产生长度为 126，输出‘0’和‘1’的 M 序列值。为保证每一个码值都能有效的输入抽气系统，程序将M序列码值以一定时间间隔逐点连续输出。

选取M序列幅值A0=±10 mm/min，在原设定运行速度20 mm/min基础上，当M序列产生码值‘1’时 赋 30 mm/min，码 值 ‘0’时 赋 10 mm/min。 在LabVIEW程序中将速度命令按“S3××W”格式以字符串形式发送至下位机，其中“××”即为赋值的运行速度值。长度为126的M序列01码值和生成程序如图6和图7所示。

图6 M序列01码幅值图

图7 M序列产生程序

2.2 直流分量处理

对采样得到的抽气系统响应信号进行频谱分析，结果表明，输出信号中包含有直流分量。为了使建立的模型能充分反映系统的动态特性，因此必须把信号中的直流分量去除。图8显示了在输入M序列后抽气系统响应信号在去除直流分量前的波形。将采集得到的输出时域信号用快速傅里叶变换（fast Fourier transformation，简称 FFT）法变换到频域信号，将信号中频率为0的分量置为0后，再将信号进行傅里叶逆变换，即可得到不含直流分量的输出信号，再经过滤波和平滑处理得到图9所示的波形。此时的输出信号值用作系统数学模型的参数辨识。

图8 响应信号去除直流分量前波形和FFT波形

图9 响应信号去除直流分量后波形和FFT波形

3 梯度校正法

梯度校正法进行参数辨识的基本思想是沿着准则函数的负梯度方向，逐渐修正模型参数的估计值，直至准则函数达到最小值[13-15]。算法结构是在旧的估计值基础上，将增益矩阵和新息的乘积求和而得到新的估计值。梯度校正法具有算法简单、计算量小等优点。

已知抽气系统的模型阶次na=nb=2，真时延d=2，在不考虑噪声的情况下，设抽气系统模型方程的差分形式为：

将式（4）改写成矩阵形式为：

其中：

式中：y(k)——抽气系统的输出值；

u(k)——输入值；

k——采样次数；

a1，a2，b1，b2——待辨识参数；

h(k)——测量矩阵；

θ——参数矩阵。

梯度校正法的目标就是根据输入输出数据来确定参数 θ和h(k)，并以此确定估计值使得准则函数J(θ)沿着负梯度方向不断修正的值，直到J(θ)达到最小值。

其中：

R(k)——N维加权矩阵，且N=4；

ΛN(k)——对角阵元素。

将去除直流分量后的输出数据按梯度校正法进行参数辨识，辨识算法流程如图10所示。使用Matlab编写模型辨识程序，程序开始执行后，程序会先给待辨识参数 θ和h(k)赋初值；按照式（8）计算加权矩阵R(k)；下一步根据J(θ)表达式和表达式计算出估计值；之后计算出模型输出值y(k)，再计算系统输出值y(k)与模型输出值ym(k)之间的偏差ey；之后再计算准则函数的梯度，判断是否达到设定循环次数，如果没有则回到权矩阵计算那一步，直到程序运行到设定的循环次数；达到设定循环次数后，将各项数值输出，绘制图形，结束程序。

图10 梯度校正法参数辨识算法流程图

待辨识参数梯度校正变化情况如图11所示，系统输出值与模型输出值偏差如图12所示。图11表明，当辨识过程进行到11次时，待辨识参数基本达到稳定状态，图12也表明此时系统输出值与模型输出值的误差基本达到稳定状态。参数辨识结果显示，待辨识参数a1=0.1185，a2=0.0569，b1= −0.013 2，b2=0.0031。由此得到抽气系统的数学模型为：

图11 待辨识参数梯度校正变化情况

图12 系统输出值与模型输出值偏差

4 数学模型的实验验证

抽气系统是燃气热值测量装置的组成部分，位于石英玻璃燃烧器的下游，由尾气收集装置、管道、测控器件等部分组成，图13是尾气收集装置和石英玻璃燃烧器的实物图。图中左侧的尾气收集装置，通过伺服电机驱动运动块上下移动的方式将燃烧室排出的尾气抽入缸内。系统工作时只需要控制尾气收集装置内活塞移动速度即可控制玻璃燃烧室进、排气口两端压差。

图13 尾气收集装置（左）和石英玻璃燃烧器（右）实物图

为验证数学模型的准确性，分别进行了两组抽气实验，第一组是抽气系统同时通入甲烷、氧气、氩气进行抽气，第二组是供气接口断开抽取空气，之后分别将两组实验中抽气时系统的输出与数学模型仿真的输出进行对比。图14是用LabVIEW中控制与仿真模块下的PID子模块结合抽气系统数学模型编写程序进行仿真。

图14 数学模型LabVIEW仿真程序

在抽气实验中，设定条件如下：

1）以系统开始供气作为0时刻，供气时使用甲烷、氧气、氩气三种气体，其体积流量比例为65∶156∶240。

2）尾气收集装置运动块运行的初始速度设为20 mm/min。

3）将期望的压差值设为35 Pa，PID控制器中设置比例系数Kp=4.2219，积分系数Ki=0.01407，微分系数Kd=0.000282。

三种气体供气抽气实验和抽空气实验输出压差变化情况与LabVIEW仿真程序输出压差变化情况分别如图15、图16所示。从两图可以看出，抽气实验输出和数学模型仿真输出均在20 s左右达到设定值附近。从抽气实验输出波形和数学模型输出波形达到较为稳定时来看，二者输出波形较为接近。三种气体供气抽气实验输出波形出现波动，是因为供气不稳导致，但输出波形仍在设定压差上下波动。模型输出和抽气实验的输出结果对比表明，数学模型能比较准确的反映系统输入输出情况。

图15 数学模型输出压差和三种气体抽气实验输出压差

图16 数学模型输出压差和抽空气实验输出压差

图17反映抽空气实验和三种气体抽气实验工况下的模型误差变化情况。抽空气时，数学模型输出结果与抽气实验输出结果的平均相对误差为0.42%，压差值平均绝对误差约为0.15 Pa；三种气体抽气实验时，数学模型输出结果与抽气实验结果的平均相对误差为0.45%，压差值平均绝对误差约为0.16 Pa。

图17 模型相对误差变化图

5 结束语

以梯度校正法为基本原理，理论和实验相结合，在LabVIEW和Matlab等工具软件的辅助下建立了燃气热值测量装置的抽气系统数学模型。以梯度校正法的数学原理设计辨识程序，具有代码量小、计算速度快、通用程度高等优点。文中首次通过建立抽气系统数学模型的方式研究燃气热值测量装置抽气控制问题；以燃烧室进、排气口压差作为反馈信号，实现了抽气系统的恒压控制。抽气系统数学模型仿真结果表明，在PID控制下模型输出结果与实际抽气控制传感器测得结果基本吻合。在抽空气实验中，数学模型与实验输出平均相对误差为0.42%；在抽三种气体实验中，数学模型与实验输出平均相对误差为0.45%。表明文中所采用梯度校正法建立抽气系统数学模型的正确性，能满足实际控制需求，为以后设计抽气系统的数字控制器提供了理论依据。