图1 分段线性化及“过冲”现象
Fig.1 Piecewise linear representation of the magnetization curve and overshooting phenomena
摘要 电磁暂态仿真中，减小仿真步长可有效抑制非线性变压器仿真过冲失真问题，但大容量变压器的励磁涌流过程往往衰减缓慢，势必会造成仿真效率低，使得仿真时间明显增加。该文建立了一种含有非线性励磁支路的大容量变压器动态相量多尺度电磁暂态模型，并提出一种动态相量磁链预测方法，通过不同仿真步长对磁链和电流进行精准预测和仿真，解决了分段线性化在大步长下出现的过冲失真问题，实现了励磁涌流宽频特征下依旧可以在不牺牲较多谐波的情况下采用大步长仿真的功能。当励磁涌流刚发生时，磁链动态相量变化较大，利用较小仿真步长对磁链进行预测，从而根据电流-磁链曲线拟合关系求得励磁涌流；当检测到磁链频谱集中在工频附近时，磁链的动态相量形式将保持恒定，采用大仿真步长预测磁链，从而根据电流-磁链关系求得准确的励磁电流，在采用大步长的同时解决了过冲失真问题，增加了仿真的计算效率。最后的仿真结果验证了该文提出的变压器模型及动态相量预测法的正确性和有效性。
关键词：动态相量预测法 励磁涌流 非线性 过冲现象
因电磁暂态仿真能有效呈现电力系统高频暂态特性，是分析电力系统动态特性的最有效的工具[1-2]。大型电力变压器具有成本高、需连续运行的特点，一般不适用于高成本的物理模拟，采用电磁暂态数字仿真呈现变压器非线性特性一直是公认有效的方法。而随着大规模交直流输电系统的发展，直流输电系统一极投入运行时，换流变压器会出现空载合闸现象，以及另一极变压器出现和应涌流，直流输电系统可能导致主网交流变压器出现严重饱和，从而引发励磁涌流，使得高频谐波馈入系统，有可能导致谐波保护误动作[3-7]。变压器铁心具有饱和特性，因此变压器电磁暂态仿真模型中含有非线性电感元件。文献[8]提出基于变压器J-A磁滞模型的高效电磁暂态模型，分析了不同仿真步长下J-A模型方程的数值特性，采用多速率仿真方法，使得仿真效率提高了1.6倍。分段线性化方法因为其简单实用的特点成为被广泛采用的非线性建模仿真方法，然而该方法在数值计算中采用较大步长容易造成过冲，从而引起失真现象。文献[9-12]分析了分段线性化变压器饱和特性，但是没有讨论如何避免采用大步长时分段线性化过冲失真问题。文献[13]建立了含有非线性励磁支路的变压器仿真模型，利用预测校正法来抑制过冲，使得仿真步长从1μs增大至50μs，但是否适用于更大的仿真步长并未讨论。文献[14-16]分析了考虑多重电力电子开关器件动作的数值振荡问题，提出基于后退欧拉法的电磁暂态仿真方法，通过计算步长的插值倒退，可以有效抑制振荡问题，但是随着倒退步骤的增多，仿真效率相比传统电磁暂态仿真并未得到较大提升。
近几年，基于移频理论的自适应步长仿真方法得到广泛关注。为加快仿真速度，文献[17-28]分别从动态相量法、移频法、自适应移频法对电力系统电磁暂态仿真进行改善，使得仿真步长可扩大至ms级，大大提高了仿真效率。移频理论和自适应移频理论由J. Marti[25-26]和K. Strunz[28]提出，其主要思想是对传统电磁暂态仿真所得实数瞬时值信号求取其解析信号，并对解析信号的频谱移动一个工频的固定频率，因此在移频后原始信号的频谱大大降低，从而根据香农采样定理可以采用大步长仿真。而这些文献中，变压器模型仅采用传统的理想变压器模型，忽略了对变压器非线性部分的建模，且其基于希尔伯特变换的移频方法因其线性变换特点无法处理非线性问题[29]，本文所用方法可以有效解决此问题。
本文提出一种变压器非线性多尺度建模方法，能够在较大时间步长仿真情况下准确刻画励磁涌流非线性特征。提出了一种基于动态相量的预测法以对发生励磁涌流时的磁链进行准确预测，既能采用毫秒级大步长仿真，又能避免分段线性化过冲失真现象，从而保证了仿真的高效性和准确性。本文研究思路为：以传统电流源的方法进行变压器非线性铁心建模，采用包含所有谐波成分的增广动态相量推导出线性电感伴随模型。电流源的求取，首先采用代表0、1、2、3次谐波的磁链动态相量对磁链进行预测，根据动态相量的慢变特性，可以采用较大步长进行准确预测，从而根据流过铁心的电流与磁链的曲线拟合关系求出励磁电流。所以本文方法可以解决励磁涌流宽频情况下，在损失较少谐波成分的同时采用大步长仿真的难题。
分段线性化作为传统处理非线性仿真的方法具有简单实用、计算量小的特点，但是如果采用较大仿真步长则容易产生过冲引起失真。磁链-电流关系满足非线性特征，非线性磁化励磁曲线如图1a所示。过冲是指当采用较大仿真步长时，磁链采用线性化求取时与真实值出现较大偏差，如图1b中的虚线所示。避免过冲问题的办法为减小仿真步长，使得c点无限接近x点，但却大大降低了仿真效率。因此，大步长情况下准确预测出磁链对避免过冲问题起到至关重要的作用。同时，采用对实数信号的梯形积分法在实现二阶精度和效率方面占有优势，但是实数信号相较于复数域的动态相量始终为快变信号，而动态相量的慢变特点则可以有效提升仿真效率。于是本文采用动态相量法进行磁链预测。
图1 分段线性化及“过冲”现象
Fig.1 Piecewise linear representation of the magnetization curve and overshooting phenomena
本文对非线性电感铁心和漏感等变压器元器件进行增广动态相量建模，将非线性部分等效为电流源并联在变压器低压侧。所以变压器模型分为线性电感建模和电流源建模。首先对线性电感（漏感）利用增广动态相量建模，其机理表述如下。
该方法所得相量包含所有谐波分量和直流分量。考虑一个实值周期函数，它用一系列表示不同频率、具有不同振幅和相位差的无穷多个正弦波的和表示。时域波形x(τ)的傅里叶级数表示[17-21]为
（1）
式中，ωs=2π/T, T为工频周期，
为第k阶随时间变化的傅里叶系数，τ∈(t-T, t]。
（2）
复数域系数被称为动态相量，能有效呈现时域函数x(t)的第k阶谐波幅值随时间变化的特征。上述方法被广泛应用于识别准稳态运行模式下的电力电子变换器仿真，被称作移频法（Shifted Frequency Modeling, SFM）[30-31]。可以清楚看到，若式（1）得到不同频率段的动态相量，就可以利用低频段信号，甚至中频和高频信号来表示原始信号。对于非线性电感和线性电感建模，需要用到式（3）、式（4）动态相量关系式。
（3）
（4）
式中，x(t)和y(t)为任意两个实数信号。然而，式（3）只描述了k阶谐波的微分变换，针对包含多次谐波的信号，采用式（1）会很大程度地增加计算量，且无法包括所有频次的谐波，从而造成失真。本文采用增广动态相量法解决此问题。
对于在一个滑动时间窗(t-T,t]内的原始信号x(t)，可以表示为以基频为基底的形式[32]，即
（5）
式中，A1和φ1分别为基频分量的幅值和相位。除了基频以外，还包括以基频为基底的所有谐波分量和直流分量，这里记作Xr(t)，即
（6）
Xr(t)与基频动态相量的和就是本文所求的增广动态相量，可以表征原始信号稳态和暂态的所有信息。
基于增广动态相量，给出变压器基本电气元件线性电感的增广动态相量伴随离散表达式，其微分方程为
（7）
因此基于式（3），线性电感的增广动态相量表达式可推导为
（8）
对式（8）利用梯形积分法离散化得
（9）
式中，Ihist为历史电流源。式（8）和式（9）对应的伴随模型如图2所示。
图2 线性电感伴随模型
Fig.2 Companion model for linear inductor
对于铁心饱和的非线性特性建模，需要准确预测每更新一个时间步长对应的磁链大小，从而根据电流-磁链关系得到对应励磁电流，建立电流源模型。相比非线性电感来表征饱和特性，电流源模型可以简化系统结构，提高仿真效率。本文拟在电流源模型的基础上进行改进。
变压器励磁涌流-磁链关系如图3所示，其中ψM表示膝点磁链，通常通过开路试验确定其轨迹。
图3 磁化特性曲线
Fig.3 Core magnetization characteristic curve
传统的方法是利用额定激励来获取数据，额外的点需要用分段线性外推技术来表征和预测，而现在普遍更倾向于采用精度更高的曲线拟合法[33-39]。所以，本文采用的曲线拟合方程[33]为
（10）
式中，LA为气隙电感；IM为膝点磁链对应的励磁电流；
为流过铁心的瞬时磁链。LA、
、IM和
为已知量。
、
、A、B、C、D分别定义为
（11）
式中，VM为膝点电压。
当考虑变压器饱和非线性时，电气信号不再是工频窄带信号，而是具有宽频特征，特别是对于励磁涌流，包含的谐波含量相较于磁链更为明显，因此对电流进行动态相量转换变为慢变信号会使得谐波次数k较大，使得计算复杂化。而磁链在变压器饱和时候的变化主要包含非周期分量，相较于励磁涌流的宽频特征，磁链（磁通）主要包含直流分量和基频分量[36]，即
（12）
式中，α为电压合闸角；θ=ωt+α；Φm为变压器的稳态磁通幅值；Φp为变压器的偏磁，Φp=Φmcosα+Φres，Φres为变压器的剩磁。
在励磁涌流期间，本文采用动态相量预测的方法，利用较大步长准确预测每个仿真步长下的磁链，原理如下：因为磁链-电流满足曲线拟合关系，由上一步仿真步长下的磁链预测出仿真步长更新一步后的磁链，则可以求出对应励磁电流。虽然磁链在变压器饱和时含有大量的直流分量和基频分量，为保证精确性，本文选用0、1、2、3次谐波的动态相量表示磁链如式（13）所示，采用动态相量预测如式（13）～式（16）和图4所示，图中历史磁链分量对应式（13）的动态相量磁链表达式。
（13）
（14）
（15）
（16）
式中，上标“～”表示预测值。
图4 磁链预测示意图
Fig.4 Flux prediction diagram
因为动态相量
为一种k次谐波慢变信号，在预测时，可以加大仿真步长，即当检测到磁链频带偏离工频较远时，式（14）、式（15）中仿真步长Δt可以从传统的电磁暂态仿真步长（通常设置为50μs）增大为500μs，当频谱中心靠近工频，则可采用2ms大步长仿真预测磁链，从而对包含0、1、2、3次谐波的磁链信号进行大步长预测。最后根据磁链-电流关系曲线拟合公式，则可以较为准确求出励磁电流值。
历史磁链动态相量求取方法说明如下：铁心磁链与其两端的电压满足
（17）
式（17）在动态相量域可表示为
（18）
式中，k=0, 1, 2, 3，分别代表磁链的直流分量、基频分量、2次谐波、3次谐波分量；ωs=2π/T，T=0.02s。
对式（18）进行梯形积分离散化得
（19）
在假定已知饱和铁心两端电压的基础上，根据式（19）可得各倍频下的历史磁链动态相量值。
传统变压器电磁暂态仿真，需要对电压电流瞬时值进行差分，因此小步长仿真成为必要条件。而基于增广动态相量的差分模型，因为采用如式（2）所示的平均化思想，所以可以实现大步长仿真，结合磁链的有效预测，亦可避免过冲现象。针对不同励磁涌流过程，本文采用不同的仿真步长，励磁涌流不同阶段仿真步长设置见表1。
表1 励磁涌流不同阶段仿真步长设置
Tab.1 Setting of time step size of diverse inrush current stage
阶段步长/ms磁链频谱中心频率fc/Hz 稳态运行250 励磁涌流初始0.5工频附近(50±10) 励磁涌流衰减期2接近工频(50±5)
首先对理想变压器建模进行说明。单相非线性变压器如图5所示，变压器一次侧相对于二次侧的电压比为a，其中不包含低压侧并联的电流源，表达式为
（20）
式中，
为流经两个端口的增广动态相量域电流，
为两端口电压相量值，
。由
，
可得节点导纳矩阵
和历史电流源
为
（21）
（22）
式中，YL为式（9）中线性电感等效导纳值。考虑饱和的变压器模型，可以在低压侧并联一个电流源。
图5 单相非线性变压器
Fig.5 Single phase saturable transformer model
该模型表达为矩阵形式有
（23）
式中，
，
为非线性电流源。电压矩阵中的各元素关系和理想变压器中相同。一次侧和二次侧增广动态相量域电流以及
为
（24）
（25）
（26）
本文采用增广动态相量法对变压器元件建模，并根据磁链频谱特性采用不同的仿真步长进行动态相量预测，整体计算流程如图6所示。
在未发生励磁涌流的稳态期间，可采用2ms大步长仿真，发生励磁涌流后分为两个阶段：初期和衰减期，分别采用较小步长500μs和大步长2ms仿真。详细流程如下。
图6 变压器多尺度仿真方法流程
Fig.6 Multi-time scale simulation flow chart of transformer
（1）算例数据录入，形成增广动态相量域下的变压器节点导纳矩阵YT如式（21）所示。
（2）采用递归离散傅里叶变换（Recursive Discrete Fourier Transform, RDFT）[37]对一个周波（即0.02s）滑窗内的线路电流频率进行持续在线检测，判别是否发生励磁涌流，通常采用在线检测二倍频谐波分量和基波的比值来判定是否发生励磁涌流[38]。若当前时刻未发生励磁涌流或者涌流消失，则执行步骤（3）；否则跳转到步骤（4）。
（3）采用大步长2ms的动态相量预测法对磁链进行准确预测，从而利用曲线拟合式（10）得出对应励磁涌流，进而跳转到步骤（6）。
（4）采用RDFT对一个周波（即0.02s）滑窗内的变压器磁链中心频率进行持续在线监测。若远离工频频率50Hz，如偏移量大于或等于10Hz，则进行步骤（5）；否跳转至步骤（3）。
（5）采用较小步长500μs的动态相量预测式（16）对磁链进行准确预测，从而利用曲线拟合式（10）得出对应励磁涌流。
（6）根据励磁涌流，获取历史电流源矩阵IhTn，从而根据式（23）求取变压器电气量。
（7）时间更新一个步长，若小于总时长，则跳转到步骤（2）；否则仿真结束。
本节对图7所示变压器进行空载合闸仿真，设置开关断开，开关在0.13s时合闸，本算例中三相短路未设置发生。电路参数见表2。仿真工具采用Matlab编程实现。变压器二次侧接一个简单电阻来模拟二次侧经输电线路和负载接地，电阻值采用PSCAD 官网模型。
图7 变压器励磁涌流电路
Fig.7 Inrush current model of transformer
表2 变压器励磁涌流电路参数
Tab.2 Parameters of transformer inrush current circuit
参数数值 E(t)/kV R1/ΩR2/ΩR3/Ω11100 L/H0.1 变压器视在功率/(MV∙A)变压器拐点电压VM(pu)变压器空心电抗LA(pu)5001.170.1 系数K1.15＜K＜1.25 IM(%)2
对该电路所有电气元件进行增广动态相量建模仿真。磁链预测采用k=0、1、2、3阶动态相量的动态相量预测法，仿真步长取500μs和2ms。励磁涌流判定发生时，采用小仿真步长500μs的动态相量预测法对磁链进行预测，当磁链频谱中心频率fc接近工频后，励磁涌流衰减速度变缓，可采用2ms大步长对励磁电流进行预测。
为了验证磁链动态相量预测法的有效性，首先对磁链的各次谐波动态相量仿真展示如下。对空载合闸系统进行0～2s的仿真如图8所示，各次谐波的动态相量在励磁涌流发生时依旧是一个慢变信号，基于此，可以说明式（14）、式（15）在励磁涌流初始时刻选择Δt=500μs，磁链衰减过程中Δt=2ms时能满足预测精度。
图8 空载合闸期间磁链各次谐波动态相量
Fig.8 Different orders of flux DPs
将预测的动态相量代入式（16），可得预测磁链如图9所示，其中0～1.8s选择采用小步长500μs预测仿真，1.8～2s选择大步长2ms进行仿真，并和1μs仿真步长波形进行了对比。图10和图11分别为图9在0～0.55s和1.8～2s的放大图。证明了所提预测算法在较大步长（500μs和2ms）仿真下的准确性，相较于传统采用1μs的电磁暂态仿真，仿真步长分别提高了500倍和2 000倍。另外，图11加入了线性外推法[37]和本文所提预测算法作对比，在2ms的仿真步长下，本文预测算法具有更高的预测精度，磁链预测精度对比见表3。
图9 仿真步长为500μs和2ms的增广动态相量磁链预测
Fig.9 Dynamic phasor prediction of flux under simulation time step 500μs, 2ms
图10 仿真步长为500μs的磁链预测放大图
Fig.10 Zoomed-in plot of prediction with 500μs simulation time step
图11 仿真步长为2ms的磁链预测放大图
Fig.11 Zoomed-in plot of prediction with 2ms time step
表3 磁链预测精度对比
Tab.3 Comparison of flux linkage prediction accuracy
预测方法步长/ms精度（%） 本文预测算法0.599.4 299.1 线性外推法279.8
为了说明本文预测算法在加大仿真步长后的精度，选取2ms大步长仿真对整个励磁涌流发生到结束进行仿真，查看本算法在大步长下的精度。对图7所示变压器进行外部故障恢复性励磁涌流仿真，设置变压器于2.5s时刻在二次侧发生三相短路，故障持续0.05s后消除。图12和图13分别为空载合闸励磁涌流和恢复性励磁涌流发生时段内，2ms步长下的预测磁链和参考磁链的比较，可以看到，如果在牺牲一定精度的情况下，即使在励磁涌流发生初始时刻，本文所提动态相量预测法依旧可以在2ms仿真步长下较为准确地预估磁链。当采用更小的仿真步长，例如500μs时，则预测得更为准确，如图10所示。同时，由图12、图13所示，预测磁链和参考磁链之间产生了较小延时，原因为：磁链被表示为0、1、2、3次谐波的动态相量形式，相比实时信号具有慢变特点。所以在稳态时，采用式（14）、式（15）对磁链动态相量进行预测。当稳态或者趋近稳态时，即磁链频谱中心距离基频10Hz以内时，磁链动态相量为平稳的直线如图8所示。采用式（14）、式（15）进行预测则不会产生时延。但是在励磁涌流发生初始时段内，即磁链频谱中心距离基频10Hz以外，对应图8中0.13～0.2s时间段内，四个相量变化较为明显，如果此时采用式（14）、式（15）进行磁链预测，Δt选为大步长2ms，则式（14）、式（15）的约等关系不再成立，所以用该方法的预测就会产生时延。根据磁链-电流关系，磁链的延时会导致励磁涌流的延时，从而对涌流仿真结果造成影响，但是因为预测磁链的幅值最大值并没有超过参考磁链，所以励磁涌流也不会发生过冲，从而能较好地对涌流最大值进行较为准确的仿真。然而如果采用本文算法流程，在励磁涌流暂态时段，即磁链频谱中心偏移量大于10Hz采用较小步长，则可以消除延时，获取更精确的励磁涌流。
图12 2ms仿真步长空载合闸时刻磁链预测示意图
Fig.12 Dynamic phasor prediction of flux under simulation time step 2ms
图13 恢复性励磁涌流时刻磁链预测示意图
Fig.13 Dynamic phasor prediction of flux when recovery inrush current occurs under 2ms time step
图14和图15分别为采用500μs和2ms的仿真步长计算得出的励磁涌流并和1μs仿真步长标准励磁涌流之间的对比。经过对比发现，误差最大约为2.5%，在允许范围之内，励磁涌流预测精度对比见表4。这大大提高了励磁涌流的仿真效率，同时也避免了励磁涌流过冲的问题。对于大型变压器，励磁涌流衰减时间往往达到几十秒[36]，采用本文所提方法可有效提高仿真效率。
图14 仿真步长为500μs的励磁涌流对比
Fig.14 Inrush current using 500μs simulation time step
图15 仿真步长为2ms的励磁涌流
Fig.15 Inrush current using 2ms simulation time step
表4 励磁涌流预测精度对比
Tab.4 Comparison of inrush current prediction accuracy
涌流发生过程步长/ms精度（%） 暂态过程0.597.5 趋近稳态299.2
本文对励磁涌流进行大步长仿真，采用不同的仿真方法引入两类不同的误差，即分段线性化引入的误差和采用电流-磁链拟合曲线的电流源法引入的误差。后者通过对磁链进行准确预测，继而通过已有的电流-磁链拟合公式对励磁涌流进行仿真，这其中对磁链在较大仿真步长下的预测成为励磁涌流仿真误差的主要来源。关于不同方法导致的仿真误差分类如图16所示。分段线性化的过冲原理在第1节做了简要介绍，其在大步长仿真下引入的误差较大，在5.3节对比了分段线性化和本文算法的过冲百分比。电流源法中的两类磁链预测法为线性外推法和动态相量预测法，实现的功能皆为对磁链进行预测，而线性外推在仿真步长较大（2ms）时预测精度相较于本文所提动态相量预测变差，如图11所示。图11中的误差是由于仿真过程中没有得到正确的磁链而产生的误差，电流-磁链曲线采用2.1节中的曲线拟合方程如式（10）所示，当磁链ψm得到准确预测，代入式（10）中，即可得到准确的励磁涌流im。电流-磁链曲线通过工程试验拟合所得，线性外推和动态相量预测所得磁链ψm都将代入式（10）得到励磁涌流im，误差仅来源于磁链是否准确预测。在稳态时，动态相量表示下的磁链为平稳的正弦信号包络线，而作为对比的线性外推法并没有基于动态相量模型，而是正弦磁链形态，采用线性外推在大步长下误差则会增大。线性外推预测磁链误差分析可参考文献[37]。在暂态时，当采用较小仿真步长，动态相量表示的包络线始终比原始正弦信号变化较慢，基于本文所提动态相量预测的方法可在较大仿真步长（500μs）下比基于快变瞬时信号的线性外推预测精度更高。
图16 励磁涌流仿真误差来源分类
Fig.16 Errors of inrush current simulation
暂态仿真中的误差采用差值方均根进行计算，即
（27）
式中，N为一个工频周期内的采样点数量；i1,j和i2,j分别为仿真励磁涌流和参考励磁涌流第j个采样值，j从1～N进行差值方均根计算。式（27）也可用于磁链误差计算。
本节针对现有过冲抑制方法——插值后退法进行介绍，并与本文所提算法在计算精度和计算耗时上进行对比。插值后退仿真过程如图17所示，其中tn1、tn2为饱和区到退饱和区的时间点。
图17 插值后退仿真过程
Fig.17 Interpolation backward method-based simulation
在一个仿真步长中，当遇到两次变压器饱和，则需要进行6步插值后退，如图17所示，从而避免发生像图1所示的过冲现象。插值后退法通过计算仿真步长的插值倒退，即图17中的第2步和第5步，则可以大幅减小或者完全抑制过冲问题。插值后退法在仿真精度方面几乎可达到100%，对0～0.4s内变压器空载合闸期间发生励磁涌流采用插值后退仿真如图18所示，与本文算法仿真对比如图19所示，可见该方法的确可实现较高精度，虽然相比插值后退有一定的误差，但是误差在接受范围内。本文算法的优势是计算速度大为提高，是插值后退的50倍以上。
图18 500μs插值后退仿真
Fig.18 Interpolation backward method using 500μs
图19 500μs步长下插值后退与本文算法对比
Fig.19 Interpolation backward method and proposed method comparison using 500μs step
然而，由图17可看出，相比本文采用的电流源建模非线性部分，该方法基于分段线性化电感伴随模型，首先引入的分段线性化电感模型在每次发生饱和后电感值发生改变，系统结构随即发生改变，需要重新计算一次节点导纳矩阵。然而当检测到计算出来的电流从饱和区过渡到非饱和区，则需要重新形成新的伴随积分公式来代替原始的电感差分伴随模型，反映到图17中则为第2步的回退过程，使得整个系统的节点导纳矩阵在发生变压器饱和时又需要重新计算一次。当非线性设备较多，则大大增加计算量。目前流行的两种励磁涌流仿真方法和本文方法在精度和耗时方面的对比见表5。参考算例为用PSCAD搭建的空载合闸励磁涌流波形，仿真步长设置为1μs以保证准确性。由表5可看出，相比插值后退法，本文方法虽然在精度上稍有下降，但是计算效率大大增加。其中
如果设置分段线性化和插值后退的步长为500μs，本文算法采用500μs、2ms的变步长仿真，对比过冲百分比，发现分段线性化在采用大步长后过冲值超过额定值的5倍左右，而插值后退和本文算法的理论过冲值约为额定值的2%左右，本文算法的过冲百分比相较于插值后退法略高2.1%和0.6%，但是计算耗时却比插值后退法快了73倍和146倍，精度也在允许范围内。随着变压器衰减时间常数增大，本文算法则更加准确。
表5 变压器励磁涌流仿真方法对比
Tab.5 Inrush current simulation comparison
励磁涌流特性算法步长/ms精度计算耗时/s（仿真4s）过冲百分比（%） 涌流大小2.87kA衰减时间常数1s分段线性化0.001步长1μs下99.9%步长大于1μs下发生严重过冲199.28（暂态下137.11，稳态下62.17）497(暂态下487，稳态下563) 插值后退0.598.8%（暂态下98.5%，稳态下99.1%）22.03（暂态下15.23，稳态下6.8）1.2(暂态下1.5，稳态下0.9) 本文算法0.5ms、2ms变步长96.7%（暂态下96%，稳态下99.9%）0.30（暂态下0.17，稳态下0.13）3.3(暂态下4，稳态下0.1) 涌流大小2.95kA衰减时间常数15s分段线性化0.001步长1μs下99.6%步长大于1μs发生严重过冲133.67432 插值后退0.598.6%16.111.4 本文算法0.5ms、2ms变步长98%0.112.0
另外，如果对比三种算法在暂态和稳态下的步长优势，表5中精度和过冲百分比分别设置了暂态和稳态情况下三种算法的对比情况。暂态过程指磁链频谱中心距离基频10Hz以上的动态，稳态过程指中心频率偏移量在10Hz以内。对比发现，分段线性化在暂态和稳态时，采用500μs仿真步长均发生严重过冲。插值后退法相比本文算法，在暂态过程中，本文算法仿真速度是插值后退法的89.59倍，但是过冲百分比高了2.5%，可以接受。稳态过程中，本文算法仿真速度是插值后退法的52.3倍，且过冲百分比相较后退法低了0.8%。综合来看，本文算法在提高仿真效率方面更为出色。
同时，表5考虑了涌流大小和衰减快慢对仿真的影响，其中衰减时间对分段线性化方法的计算耗时的影响最为明显，仿真选用了衰减时间常数为1s和15s的例子进行对比，发现三种算法的耗时均有所下降，但是本文所提算法的计算耗时依旧是最低的，比插值后退法快了73倍和146倍，涌流衰减时间越久，本文算法耗时越短，且仿真精度越高。
针对变压器饱和进行线性化中的过冲现象和仿真效率低下的问题，本文提出动态相量预测法，有效避免了过冲现象。所得结论如下：
1）与传统动态相量法相比，本文所用增广动态相量包含了所有谐波分量和非周期分量，保证了仿真精度。同时，当考虑变压器饱和非线性时，电气信号不再是一个工频窄带信号，相比于包含较多谐波量的励磁电流，磁链包含的谐波量较少，励磁涌流发生时刻，利用动态相量预测法，采用较小仿真步长（例如500μs）对磁链进行准确预测；而信号趋近于工频窄带附近时，磁链动态相量为慢变信号，则可自适应采用较大仿真步长（例如2ms）进行预测。最后可以根据磁链-电流关系求出励磁涌流，大大提高了效率，继而解决了宽频电流信号仿真依旧可以采用大步长且谐波损失较少的问题。
2）所采用的动态相量预测法相比线性外推法，能更精确有效地在励磁涌流发生时预测磁链。本方法可有效避免采用线性化求取励磁磁链的方法所带来的过冲问题。
参考文献
[1] 姚蜀军, 刘畅, 汪燕, 等. 多频段时间尺度变换电磁暂态仿真研究[J]. 中国电机工程学报, 2019, 39(24): 7199-7208.
Yao Shujun, Liu Chang, Wang Yan, et al. A research on multi-frequency band time-scale transformation for electromagnetic transients simulation[J]. Proceedings of the CSEE, 2019, 39(24): 7199-7208.
[2] 李博, 包涌泉, 彭振东, 等. 基于改进型直流真空断路器弧后暂态仿真及介质恢复特性分析[J].电工技术学报，2021, 36(8): 1752-1760.
Li Bo, Bao Yongquan, Peng Zhendong, et al. Post-arc transient simulation and dielectric recovery analysis based on improved dc vacuum circuit breaker[J]. Transactions of China Electrotechnical Society, 2021, 36(8): 1752-1760.
[3] 吴嘉琪, 李晓华, 陈忠, 等. 考虑磁滞特性变压器PSCAD /EMTDC电磁暂态仿真建模方法及励磁差异性分析[J]. 中国电机工程学报, 2017, 37(5): 1543-1551.
Wu Jiaqi, Li Xiaohua, Chen Zhong, et al. A transformer model with hysteresis characteristics for electromagnetic transients based on PSCAD /EMTDC and excitation difference analysis[J]. Proceedings of the CSEE, 2017, 37(5): 1543-1551.
[4] 李春艳, 周念成, 王强钢, 等. 基于软启动的变压器励磁涌流抑制方法[J]. 电工技术学报, 2020, 35(17): 3640-3651.
Li Chunyan, Zhou Niancheng, Wang Qianggang, et al. A method to eliminate transformer inrush currents using soft-starter-based controlled energization[J]. Transactions of China Electrotechnical Society, 2020, 35(17): 3640-3651.
[5] 王小君, 毕成杰, 金程, 等. 电气化铁路不停电过分相电磁暂态及抑制措施研究[J]. 电工技术学报, 2021, 36(1): 191-202.
Wang Xiaojun, Bi Chengjie, Jin Cheng, et al. Research on electromagnetic transient and suppression measures for passing neutral section without power interruption of electrified railway[J]. Transactions of China Electrotechnical Society, 2021, 36(1): 191-202.
[6] 任于展, 武仕朴, 汪友华, 等. 基于暂态电流差极值的电力变压器剩磁测量方法[J].电工技术学报, 2021, 37(5): 1088-1097.
Ren Yuzhan, Wu Shipu, Wang Youhua, et al. Residual flux measurement method of power transformer based on extreme value of transient current difference[J]. Transactions of China Electrotechnical Society, 2021, 2021, 37(5): 1088-1097.
[7] Fang Xupeng, Tian Yingying, Ding Xiaokang, et al. Series-type switched-inductor Z-source inverter[J]. CES Transactions on Electrical Machines and Systems, 2020, 4(1): 53-60.
[8] 洪潮, 张庆华, 张野, 等. 基于J-A磁滞模型的高效电磁暂态变压器建模[J]. 南方电网技术, 2020, 14(2): 75-83.
Hong Chao, Zhang Qinghua, Zhang Ye, et al. Efficient electromagnetic transient modeling of transformer based on J-A hysteresis model[J]. Southern Power System Technology, 2020, 14(2): 75-83.
[9] 孙洋，黄家栋. 基于非饱和区等效瞬时电感波形特征的变压器励磁涌流鉴别方法[J]. 电网技术, 2010, 34(3): 88-92.
Sun Yang，Huang Jiadong. Identification of transformer magnetizing inrush based on waveform characteristics of equivalent instantaneous inductance in non-saturation zone[J]. Power System Technology, 2010, 34(3): 88-92.
[10] 宋九渊, 符玲, 熊思宇, 等. 基于二阶泰勒系数的励磁涌流识别方法[J]. 中国电机工程学报, 2020, 40(3): 1020-1030.
Song Jiuyuan, Fu Ling, Xiong Siyu, et al. Magnetizing inrush current identification method based on second-order Taylor derivative[J]. Proceedings of the CSEE, 2020, 40(3): 1020-1030.
[11] 徐岩, 周霏霏. 基于幅值特征的变压器励磁涌流和故障电流的识别[J]. 电网技术, 2011, 35(9): 211-215.
Xu Yan, Zhou Feifei. A method to distinguish inrush current of power transformer from fault current based on amplitude characteristics[J]. Power System Technology, 2011, 35(9): 211-215.
[12] 司马文霞, 刘永来, 杨鸣, 等. 考虑铁心深度饱和的单相双绕组变压器改进π模型[J]. 中国电机工程学报, 2018, 38(24): 7131-7140.
Sima Wenxia, Liu Yonglai, Yang Ming, et al. An improved π model for single-phase two winding transformers considering deep saturation of the iron core[J]. Proceedings of the CSEE, 2018, 38(24): 7131-7140.
[13] 朱翊, 房大中, 王庆平. 考虑变压器励磁非线性的电磁暂态仿真算法[J]. 电网技术, 2012, 36(4): 116-120.
Zhu Yi, Fang Dazhong, Wang Qingping. An approach for electromagnetic transient simulation of power transformers with nonlinear exciting branch[J]. Power System Technology, 2020, 14(2): 75-83.
[14] Faruque M O, Dinavahi V, Wilsun X. Algorithms for the accounting of multiple switching events in digital simulation of power-electronic systems[J]. IEEE Transactions on Power Delivery, 2005, 20(2): 1157-1167.
[15] Kuffel P, Kent K, Irwin G. The implementation and effectiveness of linear interpolation within digital simulation[J]. International Journal of Electrical Power, 1997, 19(4): 221-227.
[16] 赵帅, 贾宏杰, 李建设, 等.一种考虑多重开关动作的变步长电磁暂态仿真算法[J]. 电工技术学报, 2016, 31(12): 177-183.
Zhao Shuai, Jia Hongjie, Li Jianshe, et al. Variable step integration method on power system transient simulation with multiple switching events[J]. Transactions of China Electrotechnical Society, 2016, 31(12): 177-183.
[17] 姚蜀军, 屈秋梦, 蔡焱蒙, 等. 基于多频段动态相量法的MMC换流器建模方法[J]. 中国电机工程学报, 2020, 40(18): 5932-5942.
Yao Shujun, Qu Qiumeng, Cai Yanmeng, et al. A research of modeling method of modular multilevel converter based on multi-frequency bands dynamic phasor[J]. Proceedings of the CSEE, 2020, 40(18): 5932-5942.
[18] 王振浩, 李博文, 成龙, 等. 多电压等级直流配电网平均动态相量建模及稳定性分析[J]. 电力系统自动化, 2021, 45(10): 50-58.
Wang Zhenhao, Li Bowen, Cheng Long, et al. Average dynamic phasor modeling and stability analysis of multi-voltage-level dc distribution network[J]. Automation of Electric Power Systems, 2021, 45(10): 50-58.
[19] 高仕林, 宋炎侃, 陈颖, 等. 电力系统移频电磁暂态仿真原理及应用综述[J]. 电力系统自动化, 2021, 45(14): 173-183.
Gao Shilin, Song Yankan, Chen Ying, et al. Overview on principle and application of shifted frequency based electromagnetic transient simulation for power system[J]. Automation of Electric Power Systems, 2021, 45(14): 173-183.
[20] 解润生, 张国荣, 高凯, 等. 适用于频变输电线的动态相量谐波分析建模与仿真[J]. 电工技术学报, 2021, 36(15): 3200-3210.
Jie Runsheng, Zhang Guorong, Gao Kai, et al. Modeling and simulation of dynamic phasor harmonic analysis for frequency-dependent transmission line[J]. Transactions of China Electrotechnical Society, 2021, 36(15): 3200-3210.
[21] Shu Dewu, Xie Xiaorong, Yan Zheng, et al. A multi-domain co-simulation method for comprehensive shifted-frequency phasor dc-grid models and EMT ac-grid models[J]. IEEE Transactions on Power Electronics, 2019, 34(11): 10557-10574.
[22] Xia Yue, Chen Ying, Song Yankan, et al. Multi-scale modeling and simulation of DFIG-based wind energy conversion system [J]. IEEE Transactions on Energy Conversion, 2020, 35(1): 560-572.
[23] Gao Feng, Strunz K. Frequency-adaptive power system modeling for multiscale simulation of transients[J]. IEEE Transactions on Power System, 2009, 24(2): 561-571.
[24] Xia Yue, Chen Ying, Song Yankan, et al. An efficient phase domain synchronous machine model with constant equivalent admittance matrix[J]. IEEE Transactions on Power Delivery, 2019, 34(3): 929-940.
[25] Zhang Peng. Shifted frequency analysis for EMTP simulation of power system dynamics[D]. Vancouver: University of British Columbia, 2009.
[26] Marti J, Dommel H, Bonatto B, et al. Shifted frequency analysis (SFA) concepts for EMTP modelling and simulation of power system dynamics[C]//2014 Power Systems Computation Conference (PSCC), Wroclaw, Poland, 2014: 1-8.
[27] Gao Feng. Frequency-adaptive modeling for multi-scale simulation of transients in power systems[D]. Washington: University of Washington, 2008.
[28] Gao Feng, Strunz, K. Modeling of constant distributed parameter transmission line for simulation of natural and envelope waveforms in power electric networks[C]//The 37th North American Power Symposium, Ames, IA, 2005: 247-252.
[29] 陈颖. 基于宽频域移频理论的大规模交直流互联电网多时间尺度暂态建模与高性能仿真[D]. 北京: 清华大学, 2018: 9-12.
[30] Shu Dewu, Dinavahi V, Xie Xiaorong, et al. Shifted frequency modeling of hybrid modular multilevel converters for simulation of MTDC grid[J] IEEE Transactions on Power Delivery, 2018, 33(3): 1288-1298.
[31] Shu Dewu, Xie Xiaorong, Jiang Qirong, et al. A novel interfacing technique for distributed hybrid simulations combining EMT and transient stability models[J]. IEEE Transactions on Power Delivery, 2018, 33(1): 130-140.
[32] Mudunkotuwa K, Filizadeh S, Annakkage U. Development of a hybrid simulator by interfacing dynamic phasors with electromagnetic transient simulation[J]. IET Generation Transmission & Distribution, 2017, 11(12): 2991-3001.
[33] EMTDC manual[Z]. Winnipeg: Manitoba HVDC Research Center, 2003.
[34] 商莹. 电磁暂态算法比较及非线性模型研究[D]. 济南: 山东大学, 2009.
[35] Chiesa N, Høidalen H. Modeling of nonlinear and hysteretic iron-core inductors in ATP[C]//European EMTP-ATP Conference, Leon, Spain, 2008: 118-123.
[36] 张建松, 何奔腾, 张雪松. 变压器衰减励磁涌流的实用计算方法[J]. 电力系统自动化, 2005, 29(12): 57-60.
Zhang Jiansong, He Benteng, Zhang Xuesong. Approach to calculate attenuated magnetizing inrush current in transformer[J]. Automation of Electric Power Systems, 2005, 29(12): 57-60.
[37] Xia Yue, Chen Ying, Ye Hua, et al. Multi-scale induction machine modeling in the dq0 domain including main flux saturation[J]. IEEE Transactions on Energy Conversion, 2019, 34(2): 652-664.
[38] 卢王允, 曹勇敢, 李苗苗. 基于小波分析的变压器励磁涌流识别判据研究[J]. 电气技术, 2014, 4(12): 8-12.
Lu Wangyun, Cao Yonggan, Li Miaomiao. Research on transformer magnetizing inrush current identification criterion based on wavelet analysis[J]. Electrical Engineering, 2014, 4(12): 8-12.
[39] 王迎迎, 袁建生. 铁心材料饱和磁化强度的偏差对求解变压器涌流的影响（英文）[J]. 电工技术学报, 2019, 34(12): 2452-2459.
Wang Yingying, Yuan Jiansheng. Influence of the error of saturation magnetization of core material on the solution of transformer inrush current[J]. Transactions of China Electrotechnical Society, 2019, 34(12): 2452-2459.
Frequency Adaptive Simulation of Transformer Inrush Current Based on the Augmented Dynamic Phasor Method
Abstract The overshoot problem can be avoided by adopting small time step in electromagnetic simulations, and the inrush current of a large capacity transformer usually decays very slowly which will certainly results in low simulation efficiency. In this paper, an augmented dynamic phasor multi-timescale electromagnetic transient model was established for large-capacity transformer with nonlinear excitation branch, and a dynamic phasor flux prediction method was utilized. The current can be accurately predicted and simulated, which solves the overshoot problem of piecewise linearization under large time steps, and large time steps can be adopted without sacrificing more current harmonic components. When the inrush current occurred, a small-time step was adopted to predict the flux linkage, then magnetizing inrush current was obtained by current-flux curve fitting function. When the flux linkage spectrum was detected concentrated near the fundamental frequency, the flux linkage waveform will maintain its sinusoidal characteristics then a large time step can be adopted without sacrifice of simulation accuracy. In this way, the accurate inrush current was simulated utilizing current-flux curve fitting function, and the overshoot problem was solved using large step size to accelerate the simulation speed. Finally, the simulation results verified the correctness of the proposed transformer model and dynamic phasor prediction method.
keywords：Dynamic phasor prediction method, inrush current, nonlinear, overshoot
DOI:10.19595/j.cnki.1000-6753.tces.211076
中图分类号：TM743
国家自然科学基金资助项目（52007133）。
收稿日期 2021-07-12
改稿日期 2022-01-14
龚 振 男，1993年生，博士研究生，研究方向为电力系统多尺度仿真计算、新型电力系统稳定与控制。E-mail：gongzhen@whu.edu.cn
刘承锡 男，1985年生，教授，博士生导师，研究方向为电力系统仿真计算、电力系统稳定与控制。E-mail：liuchengxi@whu.edu.cn（通信作者）
（编辑 赫蕾）