基于分数阶传输线模型的轨道电路暂态分析

Transient Analysis of Track Circuit Based on Fractional Order Transmission Line Model

摘要

为了准确分析暂态信号对ZPW-2000A型轨道电路的影响，考虑传输线中由集肤效应引起的频变损耗问题，建立轨道电路分数阶多导体传输线（Multi-conductor Transmission Line，MTL）模型，针对ZPW-2000A轨道电路高频损耗下暂态响应分析，提出在时域内对轨道电路接收端电压的求解方法. 基于传输线理论建立轨道电路传输线系统模型，根据得到的模型建立分数阶传输线方程并对其进行求解. 首先，在空间域上利用紧凑有限差分法（Compact Finite Difference Method，CFD）将轨道电路分数阶传输线模型的偏微分方程组离散为常微分方程组；其次，利用G-L分数阶定义将以上方程组转化为整数阶常微分方程组；最后，利用精细积分与递归卷积相结合的方法，得到传输线上每点的电压与电流响应. 在双指数信号激励下，通过与状态变量法对比验证了该方法的准确性，两种求解方法的误差在7%以内，且本文方法耗时较短. 分析了不同暂态信号输入下轨面过电压变化规律，发现信号频率越大，轨面过电压幅值越小；道床电阻越大，轨面过电压幅值越大且信号从衰减到稳定的时间越长. 本文方法可以准确、高效地分析高频损耗下ZPW-2000A型轨道电路暂态响应.

Abstract

To accurately analyze the impact of transient signals on the ZPW-2000A track circuit， this paper considers the frequency-varying loss caused by the skin effect in transmission lines and establishes a fractional-order model of the Multi-conductor Transmission Line （MTL） in the track circuit. Aiming at the transient response analysis under high-frequency loss of the ZPW-2000A track circuit， a method for solving the voltage at the receiving end of the track circuit in the time domain is proposed. Based on transmission line theory， a model of the track circuit transmission line system is established， and a fractional-order transmission line equation is formulated and solved according to the obtained model. First， the partial differential equations of the fractional-order model of the track circuit transmission line are discretized into ordinary differential equations using the Compact Finite Difference Method （CFD） in the spatial domain. Then， the G-L fractional-order definition is utilized to transform these equations into integer-order ordinary differential equations. Finally， the voltage and current responses at each point on the transmission line are obtained by combining the precise integration method with recursive convolution. Under the excitation of double-exponential signals， the accuracy of this method is verified by comparison with the state variable method， with the error between the two solutions being within 7%， and the method presented in this paper requires less time. The variation of rail overvoltage under different transient signal inputs is analyzed， revealing that the higher the signal frequency， the smaller the amplitude of rail overvoltage； the greater the ballast resistance， the larger the amplitude of rail overvoltage and the longer it takes for the signal to decay to stability.

Graphical abstract

关键词

轨道电路 / 集肤效应 / 分数阶微积分 / 多导体传输线 / 暂态分析

Key words

track circuit / skin effect / fractional calculus / multi-conductor transmission line / transient analysis

轨道电路以钢轨作为信号传输导体，是保证高速铁路安全运营的重要基础设备之一^［1］. 我国高速铁路分布距离长、跨度大，运营环境复杂，铺设在室外的钢轨与信号设备容易受到电磁暂态信号的干扰^［2］，从而影响信号设备的正常工作，威胁到列车的行车安全，因此，准确分析轨道电路的暂态响应对轨道电路传输性能的影响具有重要意义.

多导体传输线方程求解方法主要分为数值法和解析法. 数值法又分为时域法和频域法. Kunz等^［3］提出运用时域有限差分法（Finite Difference Time Domian，FDTD）求解传输线方程，为传输线的数值求解提供了理论基础. 随着传输线端接设备的复杂化，王梓丞等^［4］、陈智慧等^［5］、王川川等^［6］提出基于戴维南等效电路法、状态变量法对端接频变负载求解，并结合传输线方程求解系统电磁脉冲响应的数值算法. 针对轨道电路的研究大多采用稳态信号，Mazloom等^［7］将传输线的求解方法应用到轨道电路的暂态响应分析中. 支永健^［8］运用状态变量法分析了轨道电路信号发送端为直流输入时的响应，但电路结构复杂时列写状态方程比较困难. 乔志超等^［9］提出建立信号系统的雷击仿真模型，分析了雷击接触网时各信号设备的过电压. 祁欢等^［10］提出用精细时程积分法分析轨道电路的暂态响应，但是未考虑高频信号输入下的参数损耗.

本文考虑高频情况下轨道电路MTL设备的高频损耗，针对由轨道电路建立了分数阶多导体传输线模型. 空间上利用CFD方法将其离散成分数阶常微分方程组，用G-L分数阶定义近似代替，将分数阶微分方程组转化为整数阶微分方程组，将精细积分与递归卷积方法相结合，计算得到传输线上每点的电压与电流响应. 将计算结果与状态变量法求得的结果进行对比，误差小于7%. 为轨道电路的暂态信号传输和抗干扰分析提供了理论参考.

1 轨道电路分数阶传输线模型建立与求解

1.1 轨道电路分数阶多导体传输线建模

国内铁路普遍采用ZPW-2000系列移频轨道电路，主要包括补偿电容、发送器、接收器、电缆模拟网络、调谐单元、匹配单元及钢轨线路^［11］.其中，补偿电容的设置采用了“等间距法”，确保轨道电路入口端信号和干扰比的稳定，使轨道电路有较高传输性能；发送器用来发送高精度和高稳定性的移频信号；接收器接收信号并进行解调；电缆模拟网络分为6段，长度分别为0.5 km、0.5 km、1 km、2 km、2 km、2× 2 km，安装在室内，用于对区段轨道电路进行调节和对铁路数字信号电缆（Special Protected Telegraph，SPT）进行补偿，总补偿距离为10 km；通过设计长度为29 m的调谐区，实现两个相邻轨道电路之间的电气隔绝；匹配变压器实现轨道电路与SPT电缆之间的连接. 根据轨道电路传输特性，将钢轨等效为均匀传输线^［12］，按照标准《ZPW-2000轨道电路技术条件》（TB/T 3206—2017），设定钢轨总长度为1 km，建立轨道电路多导体传输线等效电路，如图1所示. 其中，r₁、r₂为钢轨的自电阻，l₁₁、l₂₂为钢轨的自电感，l₁₂为钢轨间的互电感，g₁₁、g₂₂为钢轨与大地间的漏电导，g₁₂为道床电阻，一般用电导表征，c₁₁、c₂₂为钢轨的自电容，c₁₂为钢轨的互电容，Z₁、Z₂表示调谐区.

基于多导体传输线理论的轨道电路传输线方程为：

\left\{\begin{array}{l}-\frac{\partial \mathit{U}(x,t)}{\partial x}=\mathit{R}\mathit{I}(x,t)+\mathit{L}\frac{\partial \mathit{U}(x,t)}{\partial t}\\ -\frac{\partial \mathit{I}(x,t)}{\partial x}=\mathit{G}\mathit{U}(x,t)+\mathit{C}\frac{\partial \mathit{U}(x,t)}{\partial t}\end{array}\right.（1）

式中： U 为电压矩阵； I 为电流矩阵； R 为传输线单位长度电阻矩阵； G 为传输线单位长度导纳矩阵； L 为传输线单位长度电感矩阵； C 为传输线单位长度电容矩阵；x、t分别为空间与时间参量.

将式（1）转换到频域，得

\left\{\begin{array}{l}-\frac{\mathrm{d}\mathit{U}(x,s)}{\mathrm{d}x}=(\mathit{R}+s\mathit{L})\mathit{I}(x,s)=\mathit{Z}(x,s)\mathit{I}(x,s)\\ -\frac{\mathrm{d}\mathit{I}(x,s)}{\mathrm{d}x}=(\mathit{G}+s\mathit{C})\mathit{U}(x,s)=\mathit{Y}(x,s)\mathit{U}(x,s)\end{array}\right.（2）

当交变电流通过导体时，导体截面上的电流密度随着靠近导体表面距离的增大而增大，这种现象称作集肤效应. 集肤效应通过减小导体的横截面积，进而增加导体的电阻来影响传输线，是传输线建模中一个需要考虑的重要参数，尤其是在高频系统中，会产生高频损耗. 高频损耗主要由导体中的集肤效应和电介质中的极化现象引起，而极化损耗对MTL影响甚微，可以忽略不计. 文献［13］通过实验得出考虑集肤效应的传输线阻抗表达式为：

Z(s)=sL+A+B\sqrt[]{s}（3）

当高频频率达到阈值时，A=0，B=K，其中K为常数.

Z{(s)}_{\mathrm{H}\mathrm{F}}=sL+K\sqrt[]{s}（4）

当低频频率达到阈值时，A=R₀，B=0.

Z{(s)}_{\mathrm{L}\mathrm{F}}=sL+R_{\mathrm{0}}（5）

当集肤深度远小于导体半径，高频频率未达到阈值时，K\sqrt[]{s}为单位面积的集肤损耗，R₀为单位长度传输线的直流电阻. 其中K=R_{\mathrm{s}}=\frac{R_{\mathrm{0}}}{\sqrt[]{\mathrm{\pi }}\sqrt[]{f_{\mathrm{0}}}}，将最终的阻抗表达式写为矩阵形式：

\mathit{Z}(s)={\mathit{R}}_{\mathrm{0}}+s\mathit{L}+{\mathit{R}}_{\mathrm{s}}\sqrt[]{s}（6）

因此高频情况下传输线模型为：

\left\{\begin{array}{l}-\frac{\mathrm{d}\mathit{U}(x,s)}{\mathrm{d}x}=({\mathit{R}}_{\mathrm{0}}+s\mathit{L}+{\mathit{R}}_{\mathrm{s}}\sqrt[]{s})\mathit{I}(x,s)\\ -\frac{\mathrm{d}\mathit{I}(x,s)}{\mathrm{d}x}=(\mathit{G}+s\mathit{C})\mathit{U}(x,s)\end{array}\right.（7）

由文献［14］可知，高频电阻R_{\mathrm{h}\mathrm{f}}随频率平方根的增大而增大，高频内电感L_{\mathrm{h}\mathrm{f}}随频率平方根的增大而减小，分别如图2和图3所示.

钢轨阻抗因集肤效应而呈现频变特性，集肤效应影响下，不同载频对应不同的钢轨阻抗值^［15］，如 表1所示.

在零初始状态下，对传输线方程运用拉普拉斯反变换，将方程从频域转化到时域的过程中，方程中的分数阶项用Caputo分数阶微积分定义下的微分形式代替，可以得到传输线的时域分数阶微分模型，如式（8）所示.

\left\{\begin{array}{l}-\frac{\partial \mathit{u}(x,t)}{\partial x}={\mathit{R}}_{\mathrm{0}}\mathit{i}(x,t)+\mathit{L}\frac{\partial \mathit{i}(x,t)}{\partial t}+{\mathit{R}}_{\mathrm{s}}\frac{{\partial }^{\mathrm{0}.\mathrm{5}}\mathit{i}(x,t)}{\partial t^{\mathrm{0.5}}}\\ -\frac{\partial \mathit{i}(x,t)}{\partial x}=\mathit{G}\mathit{u}(x,t)+\mathit{C}\frac{\partial \mathit{u}(x,t)}{\partial t}\end{array}\right.（8）

此处，分数阶微分项表示在时域中因集肤效应而产生的高频损耗. 通过在轨道电路传输线系统中引入分数阶微积分，为轨道电路基本电气参数频变现象的描述提供了新的思路.

1.2 轨道电路分数阶传输线方程求解

式（8）采用分数阶时域算法求解，对传输线偏微分方程组进行空间上的离散，使之转化为分数阶常微分方程组.传输线的分段如图4所示.

首先，对长度为L的传输线利用3阶紧凑格式差分法将传输线分为M段，每段长为\mathrm{\Delta }x=L/M，u_x表示x=(n-\mathrm{1}/\mathrm{2})\mathrm{\Delta }x(n=\mathrm{1,2},\mathrm{3},\dots ,M)处的电压；i_x表示x=n\mathrm{\Delta }x(n=\mathrm{0,1},\mathrm{2},\dots ,M)处的电流.式（9）中，每一点的值和前、后两点的取值都有关，因此，式（9）不适用于传输线的两端.对此本文采用2阶精度的紧凑差分格式^［16］对传输线的两端进行离散，首、末端离散方程分别如式（10）和式（11）所示.

\left\{\begin{array}{l}{\alpha }_{\mathrm{1}}\left(R_{\mathrm{0}}{i}_{n+\mathrm{1}}+R_{\mathrm{s}}\frac{{\mathrm{d}}^{\mathrm{0.5}}{i}_{n+\mathrm{1}}}{\mathrm{d}{t}^{\mathrm{0.5}}}+L\frac{\mathrm{d}{i}_{n+\mathrm{1}}}{\mathrm{d}t}\right)+{\alpha }_{\mathrm{2}}\left(R_{\mathrm{0}}{i}_n+R_{\mathrm{s}}\frac{{\mathrm{d}}^{\mathrm{0.5}}{i}_n}{\mathrm{d}{t}^{\mathrm{0.5}}}+L\frac{\mathrm{d}{i}_n}{\mathrm{d}t}\right)+{\alpha }_{\mathrm{1}}\left(R_{\mathrm{0}}{i}_{n-\mathrm{1}}+R_{\mathrm{s}}\frac{{\mathrm{d}}^{\mathrm{0.5}}{i}_{n-\mathrm{1}}}{\mathrm{d}{t}^{\mathrm{0.5}}}+L\frac{\mathrm{d}{i}_{n-\mathrm{1}}}{\mathrm{d}t}\right)=-\frac{u_{n+\mathrm{1}/\mathrm{2}}-u_{n-\mathrm{1}/\mathrm{2}}}{\mathrm{\Delta }x}\\ {\alpha }_{\mathrm{1}}\left(G{u}_{n+\mathrm{3}/\mathrm{2}}+C\frac{\mathrm{d}{u}_{n+\mathrm{3}/\mathrm{2}}}{\mathrm{d}t}\right)+{\alpha }_{\mathrm{2}}\left(G{u}_{n+\mathrm{1}/\mathrm{2}}+C\frac{\mathrm{d}{i}_{n+\mathrm{1}/\mathrm{2}}}{\mathrm{d}t}\right)+{\alpha }_{\mathrm{1}}\left(G{u}_{n-\mathrm{1}/\mathrm{2}}+C\frac{\mathrm{d}{u}_{n-\mathrm{1}/\mathrm{2}}}{\mathrm{d}t}\right)=-\frac{i_{n+\mathrm{1}}-i_n}{\mathrm{\Delta }x},(n=\mathrm{1,2},\cdots ,M-\mathrm{2})\end{array}\right.（9）

\left\{\begin{array}{l}-{\alpha }_{\mathrm{4}}\left(R_{\mathrm{0}}{i}_{\mathrm{0}}+R_{\mathrm{s}}\frac{{\mathrm{d}}^{\mathrm{0.5}}{i}_{\mathrm{0}}}{\mathrm{d}{t}^{\mathrm{0.5}}}+L\frac{\mathrm{d}{i}_{\mathrm{0}}}{\mathrm{d}t}\right)-{\alpha }_{\mathrm{1}}\left(R_{\mathrm{0}}{i}_{\mathrm{1}}+R_{\mathrm{s}}\frac{{\mathrm{d}}^{\mathrm{0.5}}{i}_{\mathrm{1}}}{\mathrm{d}{t}^{\mathrm{0.5}}}+L\frac{\mathrm{d}{i}_{\mathrm{1}}}{\mathrm{d}t}\right)=\frac{u_{\mathrm{1}/\mathrm{2}}-u_{\mathrm{0}}}{\mathrm{\Delta }x}\\ -{\alpha }_{\mathrm{3}}\left(G{u}_{\mathrm{1}/\mathrm{2}}+C\frac{\mathrm{d}{u}_{\mathrm{1}/\mathrm{2}}}{\mathrm{d}t}\right)-{\alpha }_{\mathrm{1}}\left(G{u}_{\mathrm{3}/\mathrm{2}}+C\frac{\mathrm{d}{u}_{\mathrm{3}/\mathrm{2}}}{\mathrm{d}t}\right)=\frac{i_{\mathrm{1}}-i_{\mathrm{0}}}{\mathrm{\Delta }x}\end{array}\right.（10）

\left\{\begin{array}{l}-{\alpha }_{\mathrm{4}}\left(R_{\mathrm{0}}{i}_M+R_{\mathrm{s}}\frac{{\mathrm{d}}^{\mathrm{0.5}}{i}_M}{\mathrm{d}{t}^{\mathrm{0.5}}}+L\frac{\mathrm{d}{i}_M}{\mathrm{d}t}\right)-{\alpha }_{\mathrm{1}}\left(R_{\mathrm{0}}{i}_{M-\mathrm{1}}+R_{\mathrm{s}}\frac{{\mathrm{d}}^{\mathrm{0.5}}{i}_{M-\mathrm{1}}}{\mathrm{d}{t}^{\mathrm{0.5}}}+L\frac{\mathrm{d}{i}_{M-\mathrm{1}}}{\mathrm{d}t}\right)=\frac{u_M-u_{M-\mathrm{1}/\mathrm{2}}}{\mathrm{\Delta }x}\\ -{\alpha }_{\mathrm{3}}\left(G{u}_{M-\mathrm{1}/\mathrm{2}}+C\frac{\mathrm{d}{u}_{M-\mathrm{1}/\mathrm{2}}}{\mathrm{d}t}\right)-{\alpha }_{\mathrm{1}}\left(G{u}_{M-\mathrm{3}/\mathrm{2}}+C\frac{\mathrm{d}{u}_{M-\mathrm{3}/\mathrm{2}}}{\mathrm{d}t}\right)=\frac{i_M-i_{M-\mathrm{1}}}{\mathrm{\Delta }x}\end{array}\right.（11）

由4阶插值公式得到{\alpha }_{\mathrm{1}}=\mathrm{1}/\mathrm{24}，{\alpha }_{\mathrm{2}}=\mathrm{1}/\mathrm{12}，{\alpha }_{\mathrm{3}}=\mathrm{23}/\mathrm{24}，{\alpha }_{\mathrm{4}}=\mathrm{11}/\mathrm{24}.

如上所示，在对轨道电路分数阶传输线方程经过空间离散后，可以得到矩阵形式的2M+1个方程为：

\mathit{P}\frac{\mathrm{d}\mathit{y}}{\mathrm{d}t}+\mathit{Q}\mathit{y}+{\mathit{Q}}_{\mathbf{s}}\frac{{\mathrm{d}}^{\mathrm{0.5}}\mathit{y}}{\mathrm{d}{t}^{\mathrm{0.5}}}+\mathit{f}=\mathrm{0}（12）

式中：\mathit{y}为u(x,t)和i(x,t)在对空间离散后形成的未知量，即\mathit{y}=[u_{\mathrm{1}/\mathrm{2}},u_{\mathrm{3}/\mathrm{2}},\dots ,u_{\mathrm{1}/\mathrm{2}},u_{M-\mathrm{1}/\mathrm{2}},i_{\mathrm{0}},i_{\mathrm{1}},\dots ,i_M{]}^{\mathrm{T}}，\mathit{f}=[\mathrm{0},\dots ,\mathrm{0},-u_{\mathrm{0}},\dots ,u_M{]}^{\mathrm{T}}，u_{\mathrm{0}}为输入信号，传输线的末端电压为u_M. [\mathit{P},\mathit{Q},{\mathit{Q}}_{\mathbf{s}}]是[R_{\mathrm{0}},R_{\mathbf{s}},L,C,G]组成的系数矩阵. 其中：

\mathit{P}=\left[\begin{array}{cccccccccc}{a}_{\mathrm{3}}& a_{\mathrm{1}}& & & & & & & & \\ a_{\mathrm{1}}& a_{\mathrm{2}}& \ddots & & & & & & & \\ & \ddots & \ddots & \ddots & & & & \mathrm{0}& & \\ & & \ddots & a_{\mathrm{2}}& a_{\mathrm{1}}& & & & & \\ & & & a_{\mathrm{1}}& a_{\mathrm{3}}& & & & & \\ & & & & & A_{\mathrm{4}}& A_{\mathrm{1}}& & & \\ & & & & & A_{\mathrm{1}}& A_{\mathrm{2}}& \ddots & & \\ & & \mathrm{0}& & & & \ddots & \ddots & \ddots & \\ & & & & & & & \ddots & A_{\mathrm{2}}& A_{\mathrm{1}}\\ & & & & & & & & A_{\mathrm{1}}& A_{\mathrm{4}}\end{array}\right]

\mathit{Q}=\left[\begin{array}{ccccccccccc}{b}_{\mathrm{3}}& b_{\mathrm{1}}& & & & -\mathrm{1}& \mathrm{1}& & & & \\ b_{\mathrm{1}}& b_{\mathrm{2}}& \ddots & \mathrm{0}& & & \ddots & \ddots & & \mathrm{0}& \\ & \ddots & \ddots & \ddots & & & & \ddots & \ddots & & \\ & \mathrm{0}& \ddots & b_{\mathrm{2}}& b_{\mathrm{1}}& & & \mathrm{0}& \ddots & \ddots & \\ & & & b_{\mathrm{1}}& b_{\mathrm{3}}& & & & & -\mathrm{1}& \mathrm{1}\\ \mathrm{1}& & & & & B_{\mathrm{4}}& B_{\mathrm{1}}& & & & \\ -\mathrm{1}& \ddots & & \mathrm{0}& & B_{\mathrm{1}}& B_{\mathrm{2}}& \ddots & & \mathrm{0}& \\ & \ddots & \ddots & & & & \ddots & \ddots & \ddots & & \\ & & \ddots & \ddots & & & & \ddots & \ddots & \ddots & \\ & \mathrm{0}& & \ddots & \mathrm{1}& & \mathrm{0}& & \ddots & B_{\mathrm{2}}& B_{\mathrm{1}}\\ & & & & -\mathrm{1}& & & & & B_{\mathrm{1}}& B_{\mathrm{4}}\end{array}\right]

{\mathit{Q}}_{\mathrm{s}}=\left[\begin{array}{ccccccccccc}& & & & & & & & & & \\ & & & & & & & & & & \\ & & \mathrm{0}& & & & & & \mathrm{0}& & \\ & & & & & & & & & & \\ & & & & & & & & & & \\ & & & & & D_{\mathrm{4}}& D_{\mathrm{1}}& & & & \\ & & & & & D_{\mathrm{1}}& D_{\mathrm{2}}& \ddots & & & \\ & & & & & & \ddots & \ddots & \ddots & & \\ & & \mathrm{0}& & & & & \ddots & \ddots & \ddots & \\ & & & & & & & & \ddots & D_{\mathrm{2}}& D_{\mathrm{1}}\\ & & & & & & & & & D_{\mathrm{1}}& D_{\mathrm{4}}\end{array}\right]

\left\{\begin{array}{l}{a}_i=C\cdot \mathrm{\Delta }x\cdot {\alpha }_i\\ A_i=L\cdot \mathrm{\Delta }x\cdot {\alpha }_i\\ b_i=G\cdot \mathrm{\Delta }x\cdot \alpha \\ B_i=R_{\mathrm{0}}\cdot \mathrm{\Delta }x\cdot {\alpha }_i\\ D_i=R_{\mathrm{s}}\cdot \mathrm{\Delta }x\cdot {\alpha }_i\end{array}\right.

由于矩阵 P 逆矩阵存在，因此，式（12）可移项变形为：

\frac{\mathrm{d}\mathit{y}}{\mathrm{d}t}={\mathit{A}}_{\mathrm{s}}\frac{{\mathrm{d}}^{\mathrm{0.5}}\mathit{y}}{\mathrm{d}{t}^{\mathrm{0.5}}}+{\mathit{A}}_{\mathrm{0}}\mathit{y}+{\mathit{f}}_{\mathrm{p}}（13）

式中：{\mathit{A}}_{\mathrm{s}}=-{\mathit{P}}^{-\mathrm{1}}{\mathit{Q}}_{\mathrm{s}}；{\mathit{A}}_{\mathrm{0}}=-{\mathit{P}}^{-\mathrm{1}}\mathit{Q}；{\mathit{f}}_{\mathrm{p}}=-{\mathit{P}}^{-\mathrm{1}}\mathit{f}.

通过CFD差分方法，轨道电路分数阶传输线方程组由偏微分方程组转化为了常微分方程组，并且其中只包含0.5阶这一个分数阶微分项. 与传统差分法相比，CFD差分方法具有更好的收敛性，并且离散所采用的点数仅为传统差分方法的1/3.

对于离散后的分数阶微分方程组，本文采用分数阶G-L定义近似代替0.5阶分数阶项，至此，分数阶微分方程组转化为了整数阶常微分方程组，利用精细积分方法求解，并结合递归卷积，可以计算得到轨道电路传输线模型上各位置处的电压与电流响应.

根据G-L定义式：

\begin{array}{l}{D}^{\alpha }y(t)=\underset{h\to \mathrm{0}}{\mathrm{l}\mathrm{i}\mathrm{m}}\frac{\mathrm{1}}{h^{\alpha }}{\displaystyle \sum _{j=\mathrm{0}}^N}{(-\mathrm{1})}^{-j}\frac{\mathrm{\Gamma }(\alpha +\mathrm{1})}{j!\mathrm{\Gamma }(\alpha -j+\mathrm{1})}\times \\ y(t-jh)\end{array}（14）

式中：\mathrm{\Gamma }(\alpha )={\int }_{\mathrm{0}}^{\infty }{e}^{-x}{x}^{\alpha -\mathrm{1}}dx=(\alpha -\mathrm{1})！为Gamma函数，，n\in R. 当时间步长h无限小，可将极限去掉，得到近似式为：

D^{\alpha }y(t)\approx \frac{\mathrm{1}}{h^{\alpha }}y(t)+K(y,t)（15）

K(y,t)=\frac{\mathrm{1}}{h^{\alpha }}{\displaystyle \sum _{j=\mathrm{1}}^N}{(-\mathrm{1})}^{-j}\frac{\mathrm{\Gamma }(\alpha +\mathrm{1})}{j!\mathrm{\Gamma }(\alpha -j+\mathrm{1})}y(t-jh)（16）

将式（15）代入式（13），得到：

\frac{\mathrm{d}y(t)}{\mathrm{d}t}=\left({\mathit{A}}_{\mathrm{0}}+\frac{\mathrm{1}}{h^{\alpha }}{\mathit{A}}_s\right)y(t)+{\mathit{A}}_{\mathrm{s}}\cdot K(y,t)+{\mathit{f}}_{\mathrm{p}}(t)（17）

对式（17）利用1阶常微分方程组的解法可以得到：

y(t)={\mathrm{e}}^{At}y(\mathrm{0})+{\int }_{\mathrm{0}}^t{\mathrm{e}}^{A(t-\tau )}[{\mathit{A}}_s\cdot K(y,\tau )+{\mathit{f}}_p(\tau )]\mathrm{d}\tau（18）

式中：\mathit{A}={\mathit{A}}_{\mathrm{0}}+\frac{\mathrm{1}}{h^{\alpha }}{\mathit{A}}_{\mathrm{s}}；y(\mathrm{0})为状态初始值.

对于矩阵指数{\mathrm{e}}^{At}，为了防止 A 矩阵的维数过高，导致计算机运算时间过长以及超出计算精度，利用精细积分方法对其进行求解^［17］：

令T={\mathrm{e}}^{\mathit{A}\tau }，m={\mathrm{2}}^N，\mathrm{\Delta }t=\tau /m，则

\begin{array}{l}{\mathrm{e}}^{\mathit{A}\mathrm{\Delta }t}\approx \mathit{I}+\mathit{A}\mathrm{\Delta }t+\frac{\mathrm{1}}{\mathrm{2}}{\left[\mathit{A}\mathrm{\Delta }t\right]}^{\mathrm{2}}+\frac{\mathrm{1}}{\mathrm{3}!}{\left[\mathit{A}\mathrm{\Delta }t\right]}^{\mathrm{3}}+\\ \frac{\mathrm{1}}{\mathrm{4}!}{\left[\mathit{A}\mathrm{\Delta }t\right]}^{\mathrm{4}}=\mathit{I}+T_{\mathrm{a}}\end{array}（19）

所以，

T={\mathrm{e}}^{\mathit{A}\tau }={\left[{\mathrm{e}}^{\mathit{A}\mathrm{\Delta }\tau }\right]}^{{\mathrm{2}}^N}={\left[\mathit{I}+{\mathit{T}}_{\mathrm{a}}\right]}^{{\mathrm{2}}^N}={\left[\mathit{I}+{\mathit{T}}_{\mathrm{a}}\right]}^{{\mathrm{2}}^{N-\mathrm{1}}}\times {\left[\mathit{I}+{\mathit{T}}_{\mathrm{a}}\right]}^{{\mathrm{2}}^{N-\mathrm{1}}}

即经过N次分解，可得：

\left[\mathit{I}+{\mathit{T}}_{\mathrm{a}}\right]\times \left[\mathit{I}+{\mathit{T}}_{\mathrm{a}}\right]=\mathit{I}+\mathrm{2}{\mathit{T}}_{\mathrm{a}}+{\mathit{T}}_{\mathrm{a}}\times {\mathit{T}}_{\mathrm{a}}

针对式（18）中的卷积项，引入递归卷积方法，将其转化为关于时间的离散形式，t=t_{n+\mathrm{1}}时，结果分为与t_n有关和具有积分形式的两项，如式（20）所示.

y(t_{n+\mathrm{1}})={\mathrm{e}}^{\mathit{A}{t}_{n+\mathrm{1}}}y(\mathrm{0})+

{\int }_{\mathrm{0}}^{t_{n+\mathrm{1}}}{\mathrm{e}}^{\mathit{A}(t_{n+\mathrm{1}}-\tau )}[{\mathit{A}}_{\mathrm{s}}\cdot K(y,\tau )+{\mathit{f}}_{\mathrm{p}}(\tau )]\mathrm{d}\tau ={\mathrm{e}}^{\mathit{A}h}y(t_n)+

{\int }_{t_n}^{t_{n+\mathrm{1}}}{\mathrm{e}}^{\mathit{A}(t_{n+\mathrm{1}}-\tau )}[{\mathit{A}}_{\mathrm{s}}\cdot K(y,\tau )+{\mathit{f}}_{\mathrm{p}}(\tau )]\mathrm{d}\tau（20）

当步长很小时，式（20）中的积分项可以用线性插值代替. 其中，迭代公式中插值系数可以取^［18］：

\left\{\begin{array}{l}{S}_{\mathrm{0}}={\mathit{A}}^{-\mathrm{1}}({\mathrm{e}}^{\mathit{A}h}-\mathit{I});S_{\mathrm{1}}={\mathit{A}}^{-\mathrm{1}}(S_{\mathrm{0}}-h\mathit{I})\\ S_{\mathrm{2}}={\mathit{A}}^{-\mathrm{1}}(\mathrm{2}{S}_{\mathrm{1}}-h^{\mathrm{2}}\mathit{I});K(y,t_n)=\frac{\mathrm{1}}{h^{\alpha }}{\displaystyle \sum _{j=\mathrm{1}}^N}{b}_{j}y(t_{n-j})\\ b_{\mathrm{1}}=-\alpha ;b_j=\left(\mathrm{1}-\frac{\mathrm{1}+\alpha }{j}\right)b_{j-\mathrm{1}}\end{array}\right.

最终结果如式（21）所示.

\begin{array}{l}y(t_{n+\mathrm{1}})={\mathrm{e}}^{\mathit{A}h}y(t_n)+\left(S_{\mathrm{0}}-\frac{\mathrm{1}}{h^{\mathrm{2}}}{S}_{\mathrm{2}}\right)[{\mathit{A}}_{\mathrm{s}}\cdot K(y,t_n)+{\mathit{f}}_{\mathrm{p}}(t_n)]+\\ \left(\frac{\mathrm{1}}{\mathrm{2}h}{S}_{\mathrm{1}}+\frac{\mathrm{1}}{\mathrm{2}{h}^{\mathrm{2}}}{S}_{\mathrm{2}}\right)\times [{\mathit{A}}_{\mathrm{s}}\cdot K(y,t_{n+\mathrm{1}})+{\mathit{f}}_{\mathrm{p}}(t_{n+\mathrm{1}})]+\\ \left(\frac{\mathrm{1}}{\mathrm{2}{h}^{\mathrm{2}}}{S}_{\mathrm{2}}-\frac{\mathrm{1}}{\mathrm{2}h}{S}_{\mathrm{1}}\right)[{\mathit{A}}_{\mathrm{s}}\cdot K(y,t_{n-\mathrm{1}})+{\mathit{f}}_{\mathrm{p}}(t_{n-\mathrm{1}})]\end{array}（21）

由式（21）可知，在知道非齐次项f跟初始条件 y（0）的情况下，即可迭代得到轨道电路在任何时刻各点的电流与电压响应.

2 方法验证

取M=200，N=20， 选取P60型钢轨、载频信号为1 700 Hz、长度为1 km的轨道区间为例，雷击波形采用常用双指数雷电电磁信号，本文选取波头为T_{\mathrm{1}}=\mathrm{1.2}\mathrm{\mu }\mathrm{s}，波长T_{\mathrm{2}}=\mathrm{50}\mathrm{\mu }\mathrm{s}的双指数波，记为1.2/50 µs 标准雷电波^［19］，幅值为30 kA为仿真基本参数.文 献［20］提出针对结构简单的传输线方程采用状态变量法求解，然而对于复杂的电路列写状态变量方程时需要借助特有的树分析.本文根据轨道电路传输线模型，求解轨道电路接收端轨面电压，并运用状态变量法对结果进行验证.

轨道电路接收端过电压计算方法对比如图5所示. 由图5可知，本文提出的时域分数阶算法与状态变量法^［21］计算结果基本一致. 计算误差在7%以内. 与状态变量法相比，本文方法无须利用矢量匹配法对传输线的频变参数和传播系数进行拟合，模型简单，避免了前期对频变参数的处理工作，具有一定的理论创新和实际应用价值.

3 影响因素对暂态响应的影响

3.1 不同频率对轨面过电压的影响

当干扰信号施加在不同载频信号的轨道上时，有必要分析轨道电路接收端过电压的变化规律.假设轨道电路发送端和接收端阻抗匹配，分别计算ZPW-2000A型轨道电路的载频信号频率取1 700 Hz、2 000 Hz、2 300 Hz和2 600 Hz时接收端的轨面过电压，得到不同频率下的接收端轨面过电压暂态响应，如图6所示.

由图6可以看出，随着信号频率从1 700 Hz增加到2 300 Hz，轨道电路接收端的电压幅值逐渐减小，并且不同频率信号趋于零的时间大约一致. 随着信号频率增大，钢轨阻抗增大，说明钢轨的集肤效应越强，信号在轨道电路上的高频损耗也就越明显，电压幅值也就越小，图6呈现的结果符合这一规律^［22］.

3.2 道床电阻对轨面过电压的影响

道床电阻作为轨道电路主要的一次参数，直接影响轨道电路传输性能. 在受到外界环境中温度和湿度的影响后，其阻抗值会随之变化，从而对轨道电路的信号传输产生干扰. 为分析不同道床电阻对轨面过电压的影响规律，选择P60型钢轨，其阻抗值是确定值，分别取道床电阻值为0.4 Ω·km、1.0 Ω·km、2.0 Ω·km、5.0 Ω·km，分析不同道床电阻对轨面电压的影响规律，如图7所示.

由图7可以看出，在同一时刻，当其他参数保持不变时，若道床电阻增加，会导致信号在传输过程中的漏泄导纳减小，进而导致漏泄电流减少，因此轨道电路接收端过电压震荡的幅值会增大；相反，若道床电阻减小，则漏泄电流增加，因而接收端轨面过电压幅值减小. 由道床电阻及轨面电压关系式可知^［23］，道床电阻越小，衰减函数越小，因此轨面过电压曲线越平缓.

4 结 论

1）本文基于传输线理论，考虑集肤效应，对 ZPW-2000A型轨道电路中因高频损耗影响而难以准确分析暂态响应的问题，建立了暂态过电压作用下的轨道电路分数阶传输线模型，将分数阶的微积分形式引入传输线系统模型，对模型采用一种时域的分数阶解法，不需要经过时频域的转换.

2）通过对ZPW-2000A型轨道电路暂态响应的分析可知，考虑高频损耗下的轨道电路接收端轨面电压的幅值小于未考虑高频损耗时得到的轨面电压幅值，符合有损传输线规律. 轨道电路载频信号频率主要影响轨道电路接收端轨面过电压的幅值，信号频率越大，轨面过电压幅值越小，但钢轨内部的频变损耗受频率的变化影响极小. 道床电阻越大，轨面过电压幅值越大且信号从衰减到稳定的时间越长.

3）将本文计算结果与状态变量法求得的波形结果进行对比，发现二者具有一致性. 与以往的多导体传输线模型相比，本文方法避免了频变参数的时域卷积处理，更有利于在时域中表示频变参数. 同时在算法方面，无须利用矢量匹配法对传输线的频变参数和传播系数进行拟合，模型简单，避免了频变参数大量的前期处理工作.分析了高频损耗、干扰信号频率、道床电阻对轨道电路接收端过电压的影响规律，为轨道电路的暂态信号传输与抗干扰分析提供了一种新的分析方法.

参考文献

[1] 张汉花，邹军，王智新，等 ．阻抗分解法及其在无砟轨道阻抗求解中的应用［J］．电工技术学报，2017，32（12）：164-170．

[2] ZHANG H H， ZOU J， WANG Z X，et al ．Impedance decomposition method and its application in calculating the impedance of ballastless track［J］．Transactions of China Electrotechnical Society，2017，32（12）：164-170．（in Chinese）

[3] 向念文，陈维江，李成榕，等 ．高铁轨道电路隔离变压器雷电暂态模型［J］．高电压技术，2016，42（5）：1594-1599．

[4] XIANG N W， CHEN W J， LI C R，et al ．Lightning transient model of isolation transformer for high-speed railway track circuit system［J］．High Voltage Engineering，2016，42（5）：1594-1599．（in Chinese）

[5] KUNZ K S， LUEBBERS R J ．The finite difference time domain method for electromagnetics［M］．Boca Raton：CRC Press，1993：10-25．

[6] 王梓丞， 郭进， 张亚东， 等. 基于FDTD接口方法的ZPW-2000轨道电路暂态分析［J］. 西南交通大学学报， 2019， 54（1）：196-201．

[7] WANG Z C， GUO J， ZHANG Y D，et al ．Transient analysis of ZPW-2000 track circuit based on FDTD interface method［J］．Journal of Southwest Jiaotong University，2019，54（1）：196-201．（in Chinese）

[8] 陈智慧， 褚庆昕 ．任意线性集总网络FDTD建模新方法［J］．西安电子科技大学学报，2008，35（2）： 262-266．

[9] CHEN Z H， CHU Q X. Novel FDTD method for modeling arbitrary linear lumped networks［J］. Journal of Xidian University， 2008， 35（2）： 262-266．（in Chinese）

[10] 王川川， 贾锐，曾勇虎，等 ．一种频变传输线系统电磁脉冲响应的数值算法［J］．北京邮电大学学报，2020，43（2）：52-58．

[11] WANG C C， JIA R， ZENG Y H，et al ．A numerical algorithm for the transient response of a frequency-dependent transmission line system excited by EMP［J］．Journal of Beijing University of Posts and Telecommunications，2020，43（2）：52-58．（in Chinese）

[12] MAZLOOM Z， THEETHAYI N， THOTTAPPILLIL R. Indirect lightning-induced voltages along a railway catenary-track multicon ductor transmission-line system with lumped components［J］. IEEE Transactions on Electromagnetic Compatibility，2011， 53（2）： 537-539．

[13] 支永健 ．弓网电弧电磁干扰传播的若干理论研究［D］．杭州：浙江大学，2013．

[14] ZHI Y J ．Some theoretical studies on pantograph-catenary arc electromagnetic interference propagation［D］. Hangzhou：Zhejiang University，2013．（in Chinese）

[15] 乔志超，向念文，阳晋，等 ．高速铁路信号系统的雷电暂态模型研究［J］．铁道学报，2018，40（9）：90-96．

[16] QIAO Z C， XIANG N W， YANG J，et al ．Research of lightning transient model for high-speed railway signal system［J］．Journal of the China Railway Society，2018，40（9）：90-96．（in Chinese）

[17] 祁欢，张友鹏 ．基于精细时程积分法的轨道电路暂态分析［J］．控制工程， 2018，25（12）：2232-2238．

[18] QI H， ZHANG Y P ．Transient analysis for track circuit based on precise time-integration［J］．Control Engineering of China，2018，25（12）：2232-2238．（in Chinese）

[19] 吴命利 ．电气化铁道牵引网的统一链式电路模型［J］．中国电机工程学报，2010，30（28）：52-58．

[20] WU M L ．Uniform chain circuit model for traction networks of electric railways［J］．Proceedings of the CSEE，2010，30（28）：52-58．（in Chinese）

[21] 龙林， 黄翰， 何月. 重载铁路牵引网雷击仿真模型比较研 究［J］. 电气化铁道， 2021， 32（3）： 55-61．

[22] LONG L， HUANG H， HE Y ．Comparative studies on lightning simulation model for traction network of heavy haul railway［J］．Electric Railway，2021，32（3）：55-61．（in Chinese）

[23] NAHMAN N， HOLT D ．Transient analysis of coaxial cables using the skin effect approximation A+B s ［J］．IEEE Transactions on Circuit Theory，1972，19（5）：443-451．

[24] 赵斌， 贾智， 王东， 等 ．电气化铁路轨道电路钢轨互阻抗计算研究［J］．铁道学报，2021，43（8）：54-61．

[25] ZHAO B， JIA Z， WANG D，et al ．Research on calculation of mutual impedance of rails in electrified railway track circuits［J］．Journal of the China Railway Society，2021，43（8）：54-61．（in Chinese）

[26] 王梓丞 ．轨道电路传输模型与故障诊断方法研究［D］．成都：西南交通大学，2019．

[27] WANG Z C ．Research on transmission model and fault diagnosis method of track circuit［D］.Chengdu：Southwest Jiaotong University，2019．（in Chinese）

[28] CANGELLARIS A C， PASHA S， PRINCE J L，et al ．A new discrete transmission line model for passive model order reduction and macromodeling of high-speed interconnections［J］．IEEE Transactions on Advanced Packaging，1999，22（3）：356-364．

[29] 张晶，尹建华，赵晓健，等 ．基于子域精细积分法的传输线故障定位研究［J］．电工技术，2018（21）：4-6．

[30] ZHANG J， YIN J H， ZHAO X J， et al. Research on transmission line fault location based on subdomain precise integration method［J］. Electric Engineering， 2018（21）： 4-6．（in Chinese）

[31] 银花，陈宁，赵尘，等 ．分数阶导数型粘弹性结构动力学方程的数值算法［J］．南京林业大学学报（自然科学版），2010，34（2）：115-118．

[32] YIN H， CHEN N， ZHAO C，et al ．A numerical algorithm for the dynamics equation of the fractional derivative viscoelasticity structure［J］．Journal of Nanjing Forestry University （Natural Sciences Edition），2010，34（2）：115-118．（in Chinese）

[33] GAO Z W， ZHANG S， HAO L，et al ．Modeling and analysis the effects of EMP on the balise system［J］．Computers，Materials & Continua，2019，58（3）：859-878．

[34] 吴振军，王丽芳，廖承林 ．分析端接频变负载的多导体传输线FDTD新方法［J］．物理学报，2009，58（9）：6146-6151．

[35] WU Z J， WANG L F， LIAO C L ．A novel FDTD method for multi-conductor transmission lines terminating in frequency-dependent loads［J］．Acta Physica Sinica，2009，58（9）：6146-6151．（in Chinese）

[36] 欧静宁. 基于多导体传输线理论的轨道电路暂态响应分 析［D］. 兰州：兰州交通大学， 2022.

[37] OU J N. Transient response analysis of track circuit based on multi-conductor transmission line theory［D］. Lanzhou：Lanzhou Jiaotong University， 2022. （in Chinese）

[38] 张友鹏， 王东， 赵斌， 等 ．计及大地影响的轨道电路钢轨自阻抗计算研究［J］．铁道科学与工程学报，2020，17（11）：2929-2937．

[39] ZHANG Y P， WANG D， ZHAO B，et al. Research on calculation of rail self impedance of track circuit considering the influence of earth［J］．Journal of Railway Science and Engineering， 2020， 17（11）： 2929-2937．（in Chinese）

[40] 阿·米·布列也夫. 轨道电路的分析与综合［M］. 孙铭甫，译. 北京：中国铁道出版社，1981： 56-73.

[41] Брылеев A M. Analysis and synthesis of track circuit ［M］. Translated by SUN M F. Beijing：China Railway Publishing House， 1981： 56-73. （in Chinese）

基金资助

国家自然科学基金资助项目(5196070019)

National Natural Science