CN105808825A - 一种适于有源配电网动态仿真的滚动式随机投影积分方法 - Google Patents

Abstract

本发明的一种适于有源配电网动态仿真的滚动式随机投影积分方法，充分考虑了有源配电系统运行过程中的随机特性和多时间尺度特征，通过将确定性投影积分算法和混合随机数值积分算法结合，并融入滚动式多随机轨迹仿真的思想，对基于随机微分‑代数方程的有源配电系统动态仿真模型进行求解，并适于在系统故障或操作条件下的随机动态仿真计算。本发明方法与传统随机数值积分算法相比，具有近似的强收敛性和弱收敛性，并可以在相同的蒙特卡罗仿真次数下，大幅提高仿真计算效率，尤其适于含多种分布式电源和负荷的有源配电系统的随机运行模拟和分析计算。

Description

一种适于有源配电网动态仿真的滚动式随机投影积分方法

技术领域

本发明涉及一种有源配电系统动态仿真方法。特别是涉及一种适于考虑系统随机特性的含多分布式电源的有源配电系统动态仿真应用的适于有源配电网动态仿真的滚动式随机投影积分方法。

背景技术

近年来，随着分布式发电(distributed generation,DG)、储能及需求侧响应等技术在配电层面的广泛应用，传统的单向无源配电系统正逐渐过渡为由于分布式电源接入而具有双向潮流的有源配电系统。有源配电系统是具备组合控制多种分布式电源、储能和负荷等能力的复杂配电系统。随着以光伏电池、小型风力发电机为代表的DG的大量接入，再加上需求侧响应技术等先进用电技术的不断发展，使得有源配电系统与传统配电系统相比，其运行过程中的随机程度大大提高。有源配电系统随机性的来源十分广泛，包括负荷的随机变化，同步发电机转子的瞬时扰动，外界环境的随机波动(风、光、温度等)和控制元件的测量误差等，这对有源配电系统的运行模拟和计算分析提出了巨大的挑战。

基于确定性方法的仿真手段目前在电力领域获得了广泛的应用，可用于分析系统运行过程中工频电气量在大扰动下(开关操作、故障等)的动态响应特性。然而，确定性方法忽略了真实物理系统中的随机特性，因此通常会得出较为保守的仿真结论。随着有源配电系统内在的随机程度不断提高，系统运行过程中的随机扰动对动态仿真结果的影响也愈发显著，使得传统确定性仿真技术在应用于具有高随机性的现代配电系统仿真时存在明显的局限性。因此，一种考虑随机特性的有源配电系统动态仿真方法对有源配电系统的规划设计和运行控制都具有至关重要的作用。

基于随机微分方程的建模与分析技术是随机系统研究的重要手段之一。目前，随机微分方程已被广泛应用于金融学中的股价和金融资产价格的随机建模，以及多个科学和工程领域中以研究在不同随机现象影响下物理系统的特性。近年来，随机微分方程在电力系统分析中的应用越来越深入，对于传统电力系统，随机微分方程多用于对负荷的随机特性进行建模，并应用于电力系统规划、仿真与小扰动分析。随着分布式电源的广泛接入，基于随机风速和光照的现代电力系统分析技术也陆续提出，如光伏发电系统和风力发电系统的随机建模等，继而为考虑随机特性的有源配电系统随机仿真奠定了基础。

有源配电系统动态仿真本质上可归结为对动力学系统时域响应的求取，分为数学建模和模型求解两部分。首先根据元件间的拓扑关系将有源配电系统各元件的特性方程构成全系统的动态仿真模型，同时考虑系统运行过程中涉及的各种随机扰动，最终形成一组随机微分-代数方程组(stochastic differential-algebraic equation,SDAE)，然后以稳态工况或潮流解为初值，求解扰动下的数值解，即逐步求得系统状态量和代数量随时间的变化曲线。

考虑随机特性的有源配电系统全系统动态仿真模型可通过一个高维随机微分-代数方程组描述，如式(1)所示。

$\left\{\begin{array}{c}\stackrel{\·}{x}=f(x,y,\mathrm{\η},\stackrel{\·}{\mathrm{\η}})\\ 0=g(x,y,\mathrm{\η})\\ \stackrel{\·}{\mathrm{\η}}=\mathrm{\α}(x,y,\mathrm{\η})+\mathrm{\β}(x,y,\mathrm{\η})\mathrm{\ξ}\end{array}\right.---\left(1\right)$

式中，为微分方程，描述有源配电系统运行的动态过程；为代数方程，表示系统潮流约束； 分别为随机微分方程的偏移量和扩散量；向量为系统状态变量，代表电机转子转速、电力电子装置的控制参数等；为代数变量，表征母线电压幅值和相角；为随机微分方程的状态变量，ξ＝dW/dt为高斯白噪声过程，其中，表示标准维纳过程。当式(1)考虑变量x，y和η的初值时，即可得到SDAE的初值形式，即随机微分方程的初值问题，如下式所示

$\left\{\begin{array}{cc}\stackrel{\·}{x}=f(x,y,\mathrm{\η},\stackrel{\·}{\mathrm{\η}})& x\left(t_0\right)=x_0\\ 0=g(x,y,\mathrm{\η})& y\left(t_0\right)=y_0\\ \stackrel{\·}{\mathrm{\η}}=\mathrm{\α}(x,y,\mathrm{\η})+\mathrm{\β}(x,y,\mathrm{\η})\mathrm{\ξ}& \mathrm{\η}\left(t_0\right)={\mathrm{\η}}_0\end{array}\right.---\left(2\right)$

此时，有源配电系统随机动态仿真可以等效为求解一个随机微分方程的初值问题。

与确定性动态仿真不同，有源配电系统随机仿真主要关注系统变量多条动态轨迹的统计特性，因此，以蒙特卡罗仿真为代表的采样方法是求解随机初值问题最常用的计算方法，并需要通过特定的数值算法对其中的每条轨迹进行求解。由于蒙特卡罗仿真的误差收敛情况与仿真次数呈正相关关系，为保证数值精度，往往需要经过大量的仿真计算。另外，有源配电系统具有较为显著的多时间尺度特征，即其仿真模型具有较强的刚性特征。因此，在有源配电系统随机动态仿真过程中，一方面需要在满足仿真精度的前提下，提高计算效率，减少蒙特卡罗仿真用时；另一方面，又需要考虑系统模型的刚性，避免出现数值稳定问题。

可见，提出一种数值精度高、计算效率高、同时考虑系统随机特性及多时间尺度特征的有源配电系统动态仿真方法，可以更加准确、高效地模拟有源配电系统在实际外部环境下的动态运行特性，对未来有源配电系统的发展和大规模广泛应用具有非常重要的意义。

发明内容

本发明所要解决的技术问题是，提供一种数值精度高、计算效率高、同时考虑系统随机特性及多时间尺度特征的适于有源配电网动态仿真的滚动式随机投影积分方法。

本发明所采用的技术方案是：一种适于有源配电网动态仿真的滚动式随机投影积分方法，包括如下步骤：

1)读取待仿真的有源配电网的拓扑连接关系、元件参数、表征各随机扰动量的随机微分方程的参数以及仿真计算参数，所述的仿真计算参数包括：仿真终止时间T，仿真步长h，蒙特卡洛仿真总次数I，随机投影算法内部积分器的积分步数k，随机投影算法外部积分器步长倍数M，并设置仿真故障及操作事件；

2)设置当前仿真次数i＝1；

3)设置滚动参数r＝0；

4)设置当前随机投影算法内部积分器的积分步数s＝1；

5)设置仿真时间t＝t+h，采用随机投影算法内部积分器计算一个步长得到有源配电系统在该时刻的状态变量x(t)和代数变量y(t)，并设置当前随机投影算法内部积分器的积分步数s＝s+1；

6)判断随机投影算法内部积分器的积分步数s是否大于k+r，是则进入下一步骤，否则返回步骤5)；

7)设置随机投影算法外部积分器的积分步长H＝Mh，设置仿真时间t＝t+H，利用随机投影算法外部积分器得到系统在该时刻的状态变量和代数变量，并设置滚动参数r＝r+1；

8)设置当前随机投影算法内部积分器的积分步数s＝1；

9)设置仿真时间t＝t+h，采用随机投影算法内部积分器计算一个步长得到有源配电系统本时刻的状态变量和代数变量，并设置当前随机投影算法内部积分器的积分步数s＝s+1；

10)根据步骤1)设置的仿真故障及操作事件，判断系统在当前时刻是否存在故障或操作，若存在，即故障或操作事件的发生时间T_event＝t，则设置滚动参数r＝r-1，并返回步骤4)，否则进入下一步骤；

11)判断随机投影算法内部积分器的积分步数s是否大于k，是则进入下一步骤，否则返回步骤9)；

12)根据步骤1)设置的仿真故障及操作事件，判断t～t+Mh时间内是否发生故障或操作，若发生，即t<T_event<t+Mh，则进入步骤14)，否则进入下一步骤；

13)设置随机投影算法外部积分器的积分步长H＝Mh，设置仿真时间t＝t+H，利用随机投影算法外部积分器得到系统在该时刻的状态变量和代数变量，转到步骤15)；

14)设置随机投影算法外部积分器的积分步长H＝T_event-t，设置仿真时间t＝T_event，利用随机投影算法外部积分器得到故障或操作发生前系统的状态变量和代数变量；

15)判断仿真时间t是否达到仿真终止时间T，若t＝T，则第i次仿真结束，设置i＝i+1，并判断此时滚动参数r是否等于k+M，若等于，则设置r＝0；若t≠T，则返回步骤8)，依据步骤8)至15)反复进行直至第i次仿真结束；

16)判断当前仿真次数是否大于蒙特卡洛仿真总次数，若大于，即i>I，则仿真结束，否则返回步骤4)，依据步骤4)至16)反复进行直至仿真结束。

所述的随机投影算法内部积分器针对有源配电系统随机模型的随机微分方程部分采用欧拉-丸山方法，对确定性部分采用基于显式改进欧拉法的微分-代数方程交替求解方法。

所述的随机投影算法外部积分器针对有源配电系统随机模型的随机微分方程部分采用高阶米尔斯坦方法，对确定性部分采用基于显式4阶龙格-库塔法的微分-代数方程交替求解方法。

本发明的一种适于有源配电网动态仿真的滚动式随机投影积分方法，充分考虑了有源配电系统运行过程中的随机特性和多时间尺度特征，基于随机微分-代数方程的有源配电系统动态仿真模型，通过确定性投影积分算法和混合随机数值积分算法的结合，并融入滚动式多随机轨迹仿真的思想，使其与传统随机数值积分算法相比，具有近似的强收敛性和弱收敛性，并可以在相同的仿真次数下，大幅提高仿真计算效率，尤其适于含多种分布式电源和负荷的有源配电系统的运行模拟和分析计算。

附图说明

图1是本发明方法的整体流程图；

图2是本发明方法的“滚动式”仿真示意图；

图3是本发明方法的强收敛性分析结果图；

图4是本发明方法的弱收敛性分析结果图；

图5是k＝3时本发明方法与传统欧拉-丸山方法弱收敛性的比较结果；

图6是低压有源配电系统算例结构图；

图中1：第一组蓄电池；2：第一组光伏电池；3：第二组光伏电池；4：第二组蓄电池；M1：中压母线；L1～L19：低压母线；Load 1～Load 8：负荷；

图7是Load 3有功功率的单轨迹仿真结果；

图8是Load 8有功功率的单轨迹仿真结果；

图9是光照强度单轨迹仿真结果；

图10是第一组光伏电池有功功率的单轨迹仿真结果；

图11是采用欧拉-丸山方法时Load 8的多轨迹仿真结果；

图12是采用随机投影积分时Load 8的多轨迹仿真结果；

图13是采用欧拉-丸山方法时第一组光伏电池有功功率的多轨迹仿真结果；

图14是采用随机投影积分时第一组光伏电池有功功率的多轨迹仿真结果；

图15是无储能且采用欧拉-丸山方法时L19母线电压多轨迹仿真结果；

图16是无储能且采用随机投影积分时L19母线电压多轨迹仿真结果；

图17是有储能且采用欧拉-丸山方法时L19母线电压多轨迹仿真结果；

图18是有储能且采用随机投影积分时L19母线电压多轨迹仿真结果。

具体实施方式

下面结合实施例和附图对本发明的一种适于有源配电网动态仿真的滚动式随机投影积分方法做出详细说明。

本发明的一种适于有源配电网动态仿真的滚动式随机投影积分方法，属于随机数值积分方法。有源配电系统接入了种类繁多的分布式电源，再加上先进用电技术的不断发展使得负荷的动态特性日趋复杂，这使得有源配电系统具有明显的随机特性。另外，有源配电系统中各元件动态响应的时间尺度差异显著，其数字仿真也显现出了较强的刚性特征，因此需要采用可同时考虑随机性和多时间尺度特征的数值积分算法以实现其仿真计算。本发明的一种适于有源配电网动态仿真的滚动式随机投影积分方法，充分考虑了有源配电系统运行过程中的随机特性和各元件的多时间尺度特征，基于随机微分-代数方程的有源配电系统动态仿真模型，通过确定性投影积分算法和混合随机数值积分算法的结合，并融入滚动式多随机轨迹仿真的思想，使其与传统随机数值积分算法相比，具有近似的强收敛性和弱收敛性，并可以在相同的仿真次数下，大幅提高仿真计算效率，尤其适于含多种分布式电源和负荷的有源配电系统运行模拟和分析计算。

本发明的一种适于有源配电网动态仿真的滚动式随机投影积分方法，是针对考虑随机特性的有源配电系统动态仿真提出的一种固定参数随机动态仿真方法，其中，求解单条仿真轨迹时随机投影算法一个投影积分步的计算步骤为：首先以步长h进行k步的小步长积分计算，随机微分方程部分采用欧拉-丸山方法，确定性部分采用基于改进欧拉法的微分-代数方程组交替求解方法；而后根据小步长仿真结果，以步长Mh进行一个大步长投影积分计算，随机微分方程部分采用高阶米尔斯坦方法，确定性部分采用基于显式4阶龙格-库塔法的微分-代数方程组交替求解方法。其中，小步长积分计算过程称为内部积分器，大步长投影积分过程称为外部积分器。而对于多轨迹蒙特卡罗仿真，则基于“滚动式”的求解思想，即将每k+M个单轨迹仿真分为一组，在每个分组内的首个投影积分步中，第j个单轨迹内部积分器的小步长积分步数为k+j-1，其后的内部积分器的积分步数恢复正常。

本发明的一种适于有源配电网动态仿真的滚动式随机投影积分方法，如图1所示，其仿真示意图如附图2所示，包括如下步骤：

1)读取待仿真的有源配电网的拓扑连接关系、元件参数、表征各随机扰动量的随机微分方程的参数以及仿真计算参数，所述的仿真计算参数包括：仿真终止时间T，仿真步长h，蒙特卡洛仿真总次数I，随机投影算法内部积分器的积分步数k，随机投影算法外部积分器步长倍数M，并设置仿真故障及操作事件；

本发明中，所述的随机投影算法内部积分器，针对有源配电系统随机模型的随机微分方程部分采用欧拉-丸山方法，即对式(2)中的随机微分方程，有

η_n+1＝η_n+hα(x_n,y_n,η_n)+β(x_n,y_n,η_n)△W_n (3)

式中，x_n，y_n和η_n分别表示t_n时刻系统的状态变量、代数变量以及随机微分方程的状态变量，ΔW_n＝W(t_n+1)-W(t_n)～N(0,h)表征标准维纳过程的随机增量。而对确定性部分则采用基于显式改进欧拉法的微分-代数方程交替求解方法。

所述的随机投影算法外部积分器针对有源配电系统随机模型的随机微分方程部分采用高阶米尔斯坦方法，即对式(2)中随机微分方程，有

$\begin{array}{c}{\mathrm{\&eta;}}_{n+k+M}={\mathrm{\&eta;}}_{n+k}+H\mathrm{\&alpha;}(x_{n+k},y_{n+k},{\mathrm{\&eta;}}_{n+k})+\mathrm{\&beta;}(x_n,y_n,{\mathrm{\&eta;}}_n)\underset{i=0}{\overset{M-1}{\&Sigma;}}{\mathrm{\&Delta;W}}_{n+k+i}\\ +\frac{1}{2}\mathrm{\&beta;}(x_n,y_n,{\mathrm{\&eta;}}_n){\mathrm{\&beta;}}^{\&prime;}(x_n,y_n,{\mathrm{\&eta;}}_n)({\left(\underset{i=0}{\overset{M-1}{\&Sigma;}}{\mathrm{\&Delta;W}}_{n+k+i}\right)}^2-H)\end{array}---\left(4\right)$

式中，β′(x_n,y_n,η_n)＝dβ(x_n,y_n,η_n)/dη_n，H为随机投影算法外部积分器的积分步长。而对确定性部分则采用基于显式4阶龙格-库塔法的微分-代数方程交替求解方法。

2)设置当前仿真次数i＝1；

3)设置滚动参数r＝0；

4)设置当前随机投影算法内部积分器的积分步数s＝1；

5)设置仿真时间t＝t+h，采用随机投影算法内部积分器计算一个步长得到有源配电系统在该时刻的状态变量x(t)和代数变量y(t)，并设置当前随机投影算法内部积分器的积分步数s＝s+1；

6)判断随机投影算法内部积分器的积分步数s是否大于k+r，是则进入下一步骤，否则返回步骤5)；

7)设置随机投影算法外部积分器的积分步长H＝Mh，设置仿真时间t＝t+H，利用随机投影算法外部积分器得到系统在该时刻的状态变量和代数变量，并设置滚动参数r＝r+1；

8)设置当前随机投影算法内部积分器的积分步数s＝1；

9)设置仿真时间t＝t+h，采用随机投影算法内部积分器计算一个步长得到有源配电系统本时刻的状态变量和代数变量，并设置当前随机投影算法内部积分器的积分步数s＝s+1；

10)根据步骤1)设置的仿真故障及操作事件，判断系统在当前时刻是否存在故障或操作，若存在，即故障或操作事件的发生时间T_event＝t，则设置滚动参数r＝r-1，并返回步骤4)，否则进入下一步骤；

11)判断随机投影算法内部积分器的积分步数s是否大于k，是则进入下一步骤，否则返回步骤9)；

12)根据步骤1)设置的仿真故障及操作事件，判断t～t+Mh时间内是否发生故障或操作，若发生，即t<T_event<t+Mh，则进入步骤14)，否则进入下一步骤；

13)设置随机投影算法外部积分器的积分步长H＝Mh，设置仿真时间t＝t+H，利用随机投影算法外部积分器得到系统在该时刻的状态变量和代数变量，转到步骤15)；

14)设置随机投影算法外部积分器的积分步长H＝T_event-t，设置仿真时间t＝T_event，利用随机投影算法外部积分器得到故障或操作发生前系统的状态变量和代数变量；

15)判断仿真时间t是否达到仿真终止时间T，若t＝T，则第i次仿真结束，设置i＝i+1，并判断此时滚动参数r是否等于k+M，若等于，则设置r＝0；若t≠T，则返回步骤8)，依据步骤8)至15)反复进行直至第i次仿真结束；

16)判断当前仿真次数是否大于蒙特卡洛仿真总次数，若大于，即i>I，则仿真结束，否则返回步骤4)，依据步骤4)至16)反复进行直至仿真结束。

下面给出具体实例：

随机微分方程的数值积分方法可通过两种收敛性描述：强收敛性和弱收敛性，其中，强收敛性表示数值方法对随机过程轨迹本身估计的准确程度。令η_n为t_n时刻由数值积分算法得到的随机微分方程准确解η(t_n)的数值解，则称算法为γ阶强收敛性的充要条件是存在一个常数C，使得

E|η_n-η(t_n)|≤Ch^γ (5)

式中，h为步长且足够小。本实例基于数值计算方法对本发明的一种适于有源配电网动态仿真的滚动式随机投影积分方法的强收敛性进行分析，并以下式所示的线性随机微分方程作为测试方程：

dη(t)＝λη(t)dt+μη(t)dW(t),η(t₀)＝η₀ (6)

式中，λ和μ为实常数；其准确解为

$\mathrm{\&eta;}\left(t\right)=\mathrm{\&eta;}\left(t_0\right)\exp \left(\right(\mathrm{\&lambda;}-\frac{1}{2}{\mathrm{\&mu;}}^2)t+\mathrm{\&mu;}W(t\left)\right)---\left(7\right)$

在测试过程中，主要关注仿真结束时刻数值解与准确解的误差，因此，令

$e_{\mathrm{\&Delta;}t}^{strong}:=E|{\mathrm{\&eta;}}_n-\mathrm{\&eta;}\left(t_n\right)|---\left(8\right)$

以表示强收敛性测试中滚动式随机投影积分算法在仿真终止时刻的误差。根据强收敛性的定义，可以得知随机投影算法为τ阶强收敛的必要条件是对于足够小的步长h，存在一个常数C，使满足

$e_{\mathrm{\&Delta;}t}^{strong}\&le;{\mathrm{C\&Delta;t}}^{\mathrm{\&tau;}}---\left(9\right)$

对式(9)两边取对数可得

${\mathrm{loge}}_{\mathrm{\&Delta;}t}^{strong}\&le;\log C+\mathrm{\&tau;}log\mathrm{\&Delta;}t---\left(10\right)$

由于可以由仿真计算得到，因此可通过最小二乘拟合求出式(10)为等式时的τ值，即该算法强收敛性阶数的估计值。

令λ＝2，μ＝1且η₀＝1，仿真时间为1s，以δt＝2^-9s为步长生成10,000个不同的离散布朗运动轨迹，对于每条轨迹分别采用5种不同的仿真步长，即h＝2^p-1δt(p取1～5的正整数)进行仿真。分别取k为1～7，M为2～10，求出相应的τ值，通过与传统随机欧拉-丸山算法比较，即可得到不同参数两种算法强收敛性的关系，如附图3所示。从图中可以看出，本发明的滚动式随机投影积分算法的强收敛性会受其参数k和M的影响，但在外部积分高阶米尔斯坦算法的作用下，在所取的算法参数区间内，随机投影算法的强收敛性整体优于传统欧拉-丸山算法。

弱收敛性则关注数值算法对随机过程统计特性估计的准确程度，由期望的误差衡量。通常，称某算法为γ阶弱收敛的充要条件是存在一个常数C，使得

|E(M(η_n))-E(M(η(t_n)))|≤Ch^γ (11)

式中，M为一个满足多项式增长条件的光滑函数。

与滚动式随机投影算法强收敛性测试类似，仍采用式(6)所示的测试方程，并令

$e_{\mathrm{\&Delta;}t}^{weak}:=|{\mathrm{E\&eta;}}_n-E\mathrm{\&eta;}\left(t_n\right)|---\left(12\right)$

表示弱收敛测试中随机投影算法在仿真终止时刻的误差。同理可得到滚动式随机投影算法与传统欧拉-丸山算法弱收敛性的关系，如附图4所示。另外，以k＝3为例，不同M值时随机投影算法与欧拉-丸山算法的弱收敛性如附图5所示。

从图中可以看出，滚动式随机投影算法的弱收敛性受参数M的影响较为明显，且随着M值的增大，其弱收敛性有所下降。与欧拉-丸山方法相比，在M值较小时，本发明的随机投影算法的弱收敛性与其近似；但当M>7时，随机投影算法的弱收敛性较欧拉-丸山方法存在较明显的差异。因此，在实际有源配电系统随机仿真应用时，当较关注仿真结果的统计特性时，为保证统计精度，建议参数M的选取不大于7。

本实例在C++编程环境中实现了本发明提出的一种适于有源配电网动态仿真的滚动式随机投影积分方法，其中，用于模拟标准维纳过程的随机增量通过C++库生成的伪随机数经过特定算法实现。并将随机投影算法与传统欧拉-丸山方法的仿真结果与计算性能进行比较以验证本发明方法的正确性和有效性。执行仿真测试的硬件平台为Intel(R)Core(TM)i7-2640CPU@2.80GHz，8GB RAM的PC机；软件环境为64位Windows 10操作系统。

本实例采用一个含分布式电源的低压有源配电系统算例(如附图6所示)对本发明方法进行测试。低压有源配电系统算例电压等级为400V，主馈线通过0.4/10kV变压器接至中压母线M1处，变压器采用常用的DYn11联结方式，低压侧设有无功补偿电容，主馈线节点间距为50m，采用三相对称线路与负荷。另外，算例中接入了多种类型的分布式电源，包括：光伏发电系统和蓄电池储能系统，各分布式电源控制方式、接入容量及有功功率输出初值如表1所示。

表1分布式电源控制方式、接入容量及有功功率初值

仿真过程中，考虑光照强度和负荷有功功率的随机变化。其中，光照强度的随机变化通过标准维纳过程描述，考虑光照强度随机性的光伏阵列模型如下式所示

$\left\{\begin{array}{c}I=N_P{I}_{ph}\left(t\right)-N_P{I}_s(e^{\frac{q}{AkT}(\frac{U}{N_S}+\frac{{\mathrm{IR}}_s}{N_P})}-1)-\frac{N_P}{R_{sh}}(\frac{U}{N_S}+\frac{{\mathrm{IR}}_s}{N_P})\\ I_{ph}\left(t\right)=\frac{S\left(t\right)}{S_{ref}}\&lsqb;I_{phref}+C_T(T-T_{ref})\&rsqb;\\ I_s=I_{sref}{\left(\frac{T}{T_{ref}}\right)}^3{e}^{\&lsqb;\frac{{\mathrm{qE}}_g}{Ak}(\frac{1}{T_{ref}}-\frac{1}{T})\&rsqb;}\\ dS\left(t\right)=s_{W}dW\left(t\right),S\left(t_0\right)=S_0\end{array}\right.---\left(13\right)$

式中，U为光伏电池输出电压；I为光伏电池输出电流；I_ph(t)为t时刻二极管扩散效应饱和电流；I_s为二极管饱和电流；q为电子电量常量，为1.602×10^-19C；k为玻尔兹曼常量，为1.381×10^-23J/K；T为光伏电池工作绝对温度值；A为二极管特性拟合系数；N_S和N_P分别为串联和并联的光伏电池数，S(t)为t时刻实际光辐照度；S_ref为标准条件下辐照度；I_phref和I_sref为标准条件辐照度下的光生电流和二极管饱和电流；T_ref为标准条件下光伏电池的工作温度；C_T为温度系数；E_g为禁带宽度；s_W为随机扰动的幅值系数，本实例中取为10。

本示例采用奥恩斯坦回归过程模拟负荷有功功率消耗的随机变化过程，其中，回归过程的随机微分方程形式如下式所示

dη(t)＝θ(μ-η(t))dt+σdW(t) (14)

式中，θ，μ和σ为参数。结合式(13)，负荷的随机动态模型如下式所示

$\left\{\begin{array}{c}{P}_L\left(t\right)=P_{L0}+{\mathrm{\&eta;}}_P\left(t\right)\\ Q_L\left(t\right)=Q_{L0}\\ {\mathrm{\&eta;}}_P\left(t\right)={\mathrm{\&theta;}}_P({\mathrm{\&mu;}}_P-{\mathrm{\&eta;}}_P(t\left)\right)dt+{\mathrm{\&sigma;}}_{P}dW\left(t\right)\end{array}\right.---\left(15\right)$

式中，P_L(t)、Q_L(t)分别表示t时刻的符合有功功率和无功功率；P_L0、Q_L0分别为负荷有功功率和无功功率初值。本实例中，各负荷有功功率初值及随机过程参数如表2所示。可以看出，Load1～Load7回归过程的均值为0，表示负荷有功功率在初值附近以方差σ_P摆动，而Load8均值为27.52kW，以模拟负荷有功功率在系统运行过程中逐渐增长的过程。

表2各负荷有功功率初值及回归过程参数

设置仿真时间为7s，分别采用确定性仿真方法，即显式4阶龙格-库塔法(RK4方法)，传统随机欧拉-丸山方法和滚动式随机投影算法方法进行仿真计算。其中，确定性RK4方法的仿真步长h为1ms；采用传统欧拉-丸山方法时，确定性部分的微分方程通过改进欧拉法求解，仿真步长h同样取1ms；采用本发明的随机投影积分方法时，参数分别取k＝3和M＝5，即内部积分步长h取1ms，外部积分步长H取5ms，维纳过程随机变量的产生步长δt与h一致。

负荷Load 3和Load 8有功功率的单轨迹仿真结果如附图7和附图8所示，光照强度和第一组光伏电池有功功率的单轨迹仿真结果如附图9和附图10所示。设置蒙特卡罗仿真次数为200次，Load 8和第一组光伏电池有功功率的随机仿真轨迹、平均值轨迹及确定性仿真结果如附图11至附图14所示。从附图9和附图10可以看出，本发明提出的滚动式随机投影算法的仿真结果与传统欧拉-丸山方法的仿真结果基本吻合，而从附图11至附图14可以看出，确定性仿真结果与随机仿真轨迹的平均值基本吻合，而采用滚动式随机投影算法时，其随机仿真轨迹和平均值轨迹与传统欧拉-丸山方法的仿真轨迹基本一致。

为进一步验证本发明提出的一种适于有源配电网动态仿真的滚动式随机投影积分方法的正确性，设置3.0s时刻L2母线发生三相经小电阻接地短路故障，短路电阻为0.5Ω，3.1s时刻故障清除，观察安装储能先后母线L19电压的变化情况，仿真结果如附图15至附图18所示。从图中可以看出，接入储能后，L19母线电压在故障期间的瞬时跌落程度具有明显改善。同样的，确定性仿真结果与随机仿真轨迹的平均值基本吻合，而采用滚动式随机投影积分算法时，其随机仿真轨迹和平均值轨迹与欧拉-丸山的仿真轨迹也基本一致。

统计附图15至附图18中故障前后电压越限的仿真轨迹个数，结果如表3所示。在不接入储能时，确定性仿真的结果不存在电压越限的情况；而在随机仿真过程中，L19母线则存在30％以上的电压越限概率。在储能接入后，电压越限的情况得到了明显改善。可以看出，确定性仿真较随机仿真的结果偏保守，在工程实际应用过程中易得出较保守的规划方案或运行策略，而随机仿真更加贴近实际情况。另外，从表中还可以看出，通过欧拉-丸山方法和随机投影积分算法得到的电压越限情况基本一致，两种方法结果的差异小于5％。

表3电压越限轨迹统计结果

传统欧拉-丸山方法和滚动式随机投影积分算法分别仿真200次的仿真用时如表4所示。在无储能的仿真场景下，欧拉-丸山方法的平均单次仿真时间约为1.475s，而随机投影积分方法的单次仿真用时约为0.734s，随机投影方法较丸山方法节约50％左右的仿真用时；同样的，在有储能的场景下，随机投影积分方法的仿真时间仍为欧拉-丸山方法的一半。因此，本发明的一种适于有源配电网动态仿真的滚动式随机投影积分方法可以在相同的仿真次数下，在保证仿真结果精度的同时，大幅提高仿真计算效率。

表4算法性能比较

综上所述，本发明提出的一种适于有源配电网动态仿真的滚动式随机投影积分方法，充分考虑了有源配电系统运行过程中的随机特性，与确定性仿真相比可以更加真实地模拟系统的实际运行状态。同时，与传统随机数值积分算法，欧拉-丸山方法相比，本发明提出的随机投影积分算法可以在相同的仿真次数下，在满足仿真结果精度的同时，大幅提高仿真速度，节约仿真用时，尤其适于含多种分布式电源和负荷的有源配电系统的运行模拟和分析计算。

Claims (3)

1.一种适于有源配电网动态仿真的滚动式随机投影积分方法，其特征在于，包括如下步骤：

1)读取待仿真的有源配电网的拓扑连接关系、元件参数、表征各随机扰动量的随机微分方程的参数以及仿真计算参数，所述的仿真计算参数包括：仿真终止时间T，仿真步长h，蒙特卡洛仿真总次数I，随机投影算法内部积分器的积分步数k，随机投影算法外部积分器步长倍数M，并设置仿真故障及操作事件；

2)设置当前仿真次数i＝1；

3)设置滚动参数r＝0；

4)设置当前随机投影算法内部积分器的积分步数s＝1；

5)设置仿真时间t＝t+h，采用随机投影算法内部积分器计算一个步长得到有源配电系统在该时刻的状态变量x(t)和代数变量y(t)，并设置当前随机投影算法内部积分器的积分步数s＝s+1；

6)判断随机投影算法内部积分器的积分步数s是否大于k+r，是则进入下一步骤，否则返回步骤5)；

7)设置随机投影算法外部积分器的积分步长H＝Mh，设置仿真时间t＝t+H，利用随机投影算法外部积分器得到系统在该时刻的状态变量和代数变量，并设置滚动参数r＝r+1；

8)设置当前随机投影算法内部积分器的积分步数s＝1；

9)设置仿真时间t＝t+h，采用随机投影算法内部积分器计算一个步长得到有源配电系统本时刻的状态变量和代数变量，并设置当前随机投影算法内部积分器的积分步数s＝s+1；

10)根据步骤1)设置的仿真故障及操作事件，判断系统在当前时刻是否存在故障或操作，若存在，即故障或操作事件的发生时间T_event＝t，则设置滚动参数r＝r-1，并返回步骤4)，否则进入下一步骤；

11)判断随机投影算法内部积分器的积分步数s是否大于k，是则进入下一步骤，否则返回步骤9)；

12)根据步骤1)设置的仿真故障及操作事件，判断t～t+Mh时间内是否发生故障或操作，若发生，即t<T_event<t+Mh，则进入步骤14)，否则进入下一步骤；

13)设置随机投影算法外部积分器的积分步长H＝Mh，设置仿真时间t＝t+H，利用随机投影算法外部积分器得到系统在该时刻的状态变量和代数变量，转到步骤15)；

14)设置随机投影算法外部积分器的积分步长H＝T_event-t，设置仿真时间t＝T_event，利用随机投影算法外部积分器得到故障或操作发生前系统的状态变量和代数变量；

15)判断仿真时间t是否达到仿真终止时间T，若t＝T，则第i次仿真结束，设置i＝i+1，并判断此时滚动参数r是否等于k+M，若等于，则设置r＝0；若t≠T，则返回步骤8)，依据步骤8)至15)反复进行直至第i次仿真结束；

16)判断当前仿真次数是否大于蒙特卡洛仿真总次数，若大于，即i>I，则仿真结束，否则返回步骤4)，依据步骤4)至16)反复进行直至仿真结束。

2.根据权利要求1所述的一种适于有源配电网动态仿真的滚动式随机投影积分方法，其特征在于，所述的随机投影算法内部积分器针对有源配电系统随机模型的随机微分方程部分采用欧拉-丸山方法，对确定性部分采用基于显式改进欧拉法的微分-代数方程交替求解方法。

3.根据权利要求1所述的一种适于有源配电网动态仿真的滚动式随机投影积分方法，其特征在于，所述的随机投影算法外部积分器针对有源配电系统随机模型的随机微分方程部分采用高阶米尔斯坦方法，对确定性部分采用基于显式4阶龙格-库塔法的微分-代数方程交替求解方法。