[发明专利]一种高效变步长数值方法

说明书

技术领域

本发明涉及一种高效的变步长数值求解方法，属于数值分析领域，可以应用于电力系统电磁暂态仿真以及计算机辅助电路分析技术等领域。

背景技术

梯形公式使用梯形近似函数的积分面积，是数值求解微分方程的常用方法；它具有二阶精度，并且是一种A稳定的方法，因此在电力系统电磁暂态仿真、计算机辅助电路分析技术等领域获得了广泛的应用。但是梯形方法也存在一些不足之处，这主要体现在：

梯形法虽然是A稳定的，但是却不是L稳定的，有可能导致数值震荡现象。例如将梯形法用于电磁暂态仿真程序，当系统中有快速的暂态过程的时候，如果仿真步长相对于暂态过程太大，非状态变量会产生围绕准确解来回震荡的现象。目前解决电磁暂态仿真程序数值震荡问题主要有两种方法：文献“Neville Watson and Jos Arrillaga,“Power System Electromagnetic Transients Simulation”Published by Institute of Engineering and Technology,London,United Kingdom,2003.”的NOS方法(Numerical Oscillation Suppression)以及文献“Marti,Jose R.；Lin,Jiming,"Suppression of numerical oscillations in the EMTP,"in Power Systems,IEEE Transactions on,vol.4,no.2,pp.739-747,May 1989.”中的CDA方法(Critical Damping Adjustment)。NOS的原理简单直观，既然变量是围绕准确的值来回震荡，那么相邻的两个计算点的中间位置必然很接近准确值，因此NOS使用半个步长的插值来抑制数值震荡。由于为了准确模拟事件发生时刻，各元件已经实现插值功能，因此通过半步长插值抑制数值震荡，比较容易实现，电磁暂态仿真程序PSCAD/EMTDC就是使用的这种方式。CDA使用半步长的后退欧拉法抑制数值震荡，由于后退欧拉法是L稳定的方法，它可以有效地抑制由于快暂态过程造成的数值震荡，而且使用半步长的后退欧拉法所得到的诺顿等值的导纳与整步长的梯形法所得到的诺顿等值的导纳是一样的，这样可以避免修改和重新分解导纳矩阵，有些电磁暂态软件(如：EMTP-RV)就是使用的这类方法。

梯形法是隐式的方法，需要在每个时步求解一个线性方程组，而且由于通常步长会影响线性方程组的系数矩阵中的元素值；所以当步长改变时，需要重新对线性方程组的系数 矩阵进行LU分解，导致计算速度的下降。在某些应用中需要却改变计算步长，例如在电磁暂态仿真程序中仿真电力电子设备时，电磁暂态仿真程序通常采用定步长仿真，而电力电子器件的动作时刻却可能落在非整数步长上，而对于电力电子器件的动作时刻的精确模拟很多情况下却是必要的。采用变步长可以方便地在电力电子器件的动作时刻插入仿真点，但是改变仿真步长会导致元件的伴随模型的导纳发生变化，这就需要重新形成系统方程的导纳阵，并重新进行LU分解，降低仿真的效率。为了解决这个问题，目前主流的电磁暂态仿真都是通过采用插值的方法在步长间插入仿真点，电磁暂态仿真程序依然使用固定的步长仿真，但是当发现某个器件在步长之间发生了事件时，通过线性插值得到事件发生时刻系统的状态。

本发明以梯形法的离散化公式为基础，并假定变量在一个步长内线性变化，通过对变量的线性插值，构建了一个新的数值方法。新数值方法可以在不改变线性方程组系数矩阵的条件下实现步长的改变，避免了线性方程组求解过程中的LU分解，大大提高了计算效率；而且新数值方法不但是A稳定的，还可以通过选取合适的参数有效地抑制计算过程中的数值震荡现象。

发明内容

为克服梯形法的不足，本发明通过对梯形离散化公式进行线性插值，得到一种新的适合于变步长的高效数值方法，新方法的描述如下：

对于如式(1)的微分方程：

梯形法的离散化公式为：

(2)式被用来从t时刻计算t+Δt时刻系统的状态，对于t时刻到t+Δt时刻间的任何时刻t′∈(t,t+Δt]，设：

t′＝t+kΔt,k∈(0,1] (3)

假设变量x(t)在t时刻到t+Δt时刻是线性变化的，则有：

将(4)和(5)式代入(2)式，整理后可得：

由(3)式可得：

Δt′＝t′-t＝kΔt (7)