CN103730900B - 电力系统多时间尺度的省地县一体化无功优化方法 - Google Patents

Abstract

本发明公开了一种电力系统多时间尺度的省地县一体化无功优化方法，本发明的多时间尺度的省地县一体化无功优化方法是一种考虑了年规律、季规律、月规律、日规律、超短时固有惯性规律和随机成份的省地县一体化无功优化方法。本发明具有能够实现多时间尺度的省地县一体化无功优化，给出多时间尺度连续的无功电压优化配置策略，实现电力系统长、中、短期各电压等级间的无功电压优化运行，提高全网的电压合格率，保障电网的持续的安全稳定运行，为制定电网的运行方式、合理且优化的利用发电资源提供理论依据，能够有效的确保电网安全稳定持续运行的特点。

Description

电力系统多时间尺度的省地县一体化无功优化方法

技术领域

本发明涉及电力系统无功分析和优化领域，尤其是涉及一种高效的电力系统多时间尺度的省地县一体化无功优化方法。

背景技术

近十几年来，电力系统安全稳定问题日益成为业界关注的焦点。在现代电网中，由于弱系统和长线路的负荷加重，电网的传输容量越来越接近其极限，使系统的稳定问题更加突出。优化电力系统的运行方式，可大幅度地改善电力系统的安全稳定环境，而电力系统的无功优化是其中的重要环节，它涉及大型、非线性、动态和离散的组合最优化问题。现有的流程是运行人员根据负荷数据、网架结构和经验设定初值，用BPA的潮流计算程序进行校验；然后人工分析BPA的潮流图，由各点电压和潮流方向人为确定各个相关的调整量。所以，即便是非常有经验的运行人员，也只能掌握住大框架，而很多细节问题不能面面俱到，要找到最优的方案是非常困难的事情。随着电力系统分析方法和计算技术的发展，急需解决上述的问题，至今未见有对传统无功优化的多方面不足进行全面解决的方案。

中国专利授权公开号：CN103199544A，授权公开日2013年7月10日，公开了一种电力系统的无功优化方法，包括以下步骤：采用牛顿-拉夫逊潮流计算程序为优化计算提供各状态变量的值和网损值；建立细菌觅食优化算法与粒子群算法结合的所述电力系统无功优化模型，确定优化目标函数；初始化细菌群体，初始化参数包括：细菌种群大小、控制变量个数、细菌位置、趋化算子次数、繁殖算子次数、迁徙算子次数、执行繁殖的比例、执行迁徙算子的概率；对细菌个体中的连续变量和离散变量进行初始化，连续变量为发电机的机端电压，离散变量包括：变压器档位、补偿电容档位；变压器档位对应变压器变比，补偿电容档位对应补偿电容量。该发明的不足之处是没有考虑时间尺度以及各区域之间无功优化局限性。

发明内容

本实用新型的发明目的是为了克服现有技术中没有根据电网电压等级不同进行不同方式无功优化以及没有考虑到不同时间尺度电网的特性进行不同无功优化方式的不足，提供了一种高效的电力系统多时间尺度的省地县一体化无功优化方法。

为了实现上述目的，本实用新型采用以下技术方案：

一种电力系统多时间尺度的省地县一体化无功优化方法，包括下列步骤：

(1-1)选择要优化的电网范围，获得对应的其电网的网架数据和电网的运行数据，计算机中设有若干个时间尺度，用户控制计算机选择至少一个1个时间尺度的网架数据和运行数据进行省地县潮流计算：

(1-1-1)读取网架数据和运行数据，网架数据包括电网网架参数，运行数据包括系统运行的参数，设定V_i＝1和δ_i＝0，设定迭代次数k＝0，最大迭代次数k_max＝M，根据网架参数生成导纳矩阵Y_B；导纳矩阵Y_B中包含电网的网架结构关系；

$Y_B=\left[\begin{array}{ccccc}{Y}_{11}& Y_{12}& Y_{13}& \mathrm{...}& Y_{1n}\\ Y_{21}& Y_{22}& Y_{23}& \mathrm{...}& Y_{2n}\\ Y_{31}& Y_{32}& Y_{33}& \mathrm{...}& Y_{3n}\\ .& .& .& .& .\\ .& .& .& .& .\\ .& .& .& .& .\\ Y_{n1}& Y_{n2}& Y_{n3}& \mathrm{..}& Y_{nn}\end{array}\right]$

式中，节点导纳矩阵的对角元Y_ii(i＝1，2，…，n)称自导纳；节点导纳矩阵的非对角元Y_ji(i＝1，2，…，n；j≠i)称互导纳。

(1-1-2)利用潮流方程 $\left\{\begin{array}{c}{P}_{Gi}-P_{Di}-\underset{j\∈S_N}{\Σ}{V}_i{Y}_{ij}{V}_j{\mathrm{cos\δ}}_{ij}=0\\ Q_{Gi}-Q_{Di}-\underset{j\∈S_N}{\Σ}{V}_i{Y}_{ij}{V}_j{\mathrm{sin\δ}}_{ij}=0\end{array}\right.,i\∈S_N,$ 计算功率不平衡量ΔP_i和ΔQ_i：

其中：i为电网中的任一个节点，j为电网中节点i之外的任一个节点，P_Gi：为节点i的发电机有功功率；Q_Gi：节点i的发电机无功功率；P_Di：节点i的负荷有功功率；Q_Di为节点i的负荷无功功率；V_i为节点i的电压幅值；Y_B是电网各个节点之间的导纳元构成的导纳矩阵，反应了电网的电气结构关系；δ_ij＝δ_i-δ_j-α_ij；δ_i为节点i的电压的相角；δ_j为节点j的电压相角；α_ij为对应导纳矩阵元素Y_ij的夹角，表示节点i和节点j之间的导纳夹角；S_N为电网节点集合；M为最大迭代次数；

(1-1-3)校验功率不平衡量ΔP_i和ΔQ_i的值，若max|(ΔP_i，ΔQ_i)|＜ε，则收敛，然后进行网络等值，其中收敛精度ε＝10^-6；

若max|(ΔP_i，ΔQ_i)|≥ε，利用 $\left[\begin{array}{cc}\frac{\∂\mathrm{\Δ}P}{\∂\mathrm{\δ}}& \frac{\∂\mathrm{\Δ}P}{\∂V}\\ \frac{\∂\mathrm{\Δ}Q}{\∂\mathrm{\δ}}& \frac{\∂\mathrm{\Δ}Q}{\∂V}\end{array}\right]\left[\begin{array}{c}\mathrm{\Δ}\mathrm{\δ}\\ \mathrm{\Δ}V\end{array}\right]=-\left[\begin{array}{c}\mathrm{\Δ}P\\ \mathrm{\Δ}Q\end{array}\right]$ 计算不平衡量ΔV和Δδ；

其中，为有功功率偏量对节点电压相角的偏导；为有功功率偏量对电压的偏导；为无功功率偏量对节点电压相角的偏导；为无功功率偏量对电压的偏导；Δδ为节点电压相角增量；ΔV为节点电压增量，ΔP为一维列向量：ΔP＝(ΔP₁，ΔP₂，…ΔP_i)，i∈S_N；ΔQ为一维列向量：ΔQ＝(ΔQ₁，ΔQ₂，…ΔQ_i)，i∈S_N；

(1-1-4)利用方程 $\left\{\begin{array}{c}\mathrm{\Δ}V=V^{(k+1)}-V^{\left(k\right)}\\ \mathrm{\Δ}\mathrm{\δ}={\mathrm{\δ}}^{(k+1)}-{\mathrm{\δ}}^{\left(k\right)}\end{array}\right.$ 求解V^(k+1)和δ^(k+1)，使k值增加1，当max|(ΔP_i，ΔQ_i)|≥ε时，重复步骤(1-1-2)；

其中，V^(k+1)为迭代次数k增加1后电网所有节点新的电压迭代值，为δ^(k+1)为迭代次数k增加1后电网所有节点新的相角迭代值，为 ${\mathrm{\δ}}^{(k+1)}=({\mathrm{\δ}}_1^{(k+1)},{\mathrm{\δ}}_2^{(k+1)},\mathrm{...}{\mathrm{\δ}}_i^{(k+1)}),i\∈S_N;$

(1-2)若k＞M，则通过最优潮流对网架参数、发电机的有功功率P_Gi、无功功率Q_Gi出力大小进行调整计算：

(1-2-1)重新读取网架数据和运行数据，设定迭代次数k₁＝0，最大迭代次数限制为k_1max＝N，设置节点电压V_i＝1和节点相角δ_i＝0，设置松弛变量l_i＝1，u_i＝1，拉格朗日乘子z_i＝1，w_i＝-0.5，y_i＝0；σ∈(0，1)称为中心参数，根据网架参数生成导纳矩阵Y_B；

(1-2-2)

当k₁＜N时，计算互补间隙Gap＝l^Tz-u^Tw，如果互补间隙Gap＜ε，则最优潮流收敛，然后进行网络等值，其中收敛精度ε＝10^-6，r为常数；

其中，l为一维列向量，l^T是l的转置，l＝(l₁，l₂，…l_i)，i∈r；u为一维列向量，u^T是u的转置，u＝(u₁，u₂，…u_i)，i∈r；z为一维列向量：z＝(z₁，z₂，…z_i)，i∈r；w为一维列向量：w＝(w₁，w₂，…w_i)，i∈r；

若Gap≥ε，则利用μ＝σ(l^Tz-u^Tw)/(2r)计算扰动因子μ，并计算修正方程组 $\{\begin{array}{c}-\[{\▿}_x^{2}f\left(x\right)-{\▿}_x^{2}h\left(x\right)y-{\▿}_x^{2}g\left(x\right)(z+w)\]\mathrm{\Δ}x+{\▿}_{x}h\left(x\right)\mathrm{\Δ}y+{\▿}_{x}g\left(x\right)(\mathrm{\Δ}z+\mathrm{\Δ}w)=L_x\\ {\▿}_x{h}^T\left(x\right)\mathrm{\Δ}x=-L_y\\ {\▿}_x{g}^T\left(x\right)\mathrm{\Δ}x-\mathrm{\Δ}x=-L_z\\ {\▿}_x{g}^T\left(x\right)\mathrm{\Δ}x+\mathrm{\Δ}u=-L_w\\ Z\mathrm{\Δ}l+L\mathrm{\Δ}z=-L_l^{\mathrm{\μ}}\\ W\mathrm{\Δ}u+U\mathrm{\Δ}w=-L_u^{\mathrm{\μ}}\end{array},$

$\left\{\begin{array}{c}{L}_x\≡{\▿}_{x}f\left(x\right)-{\▿}_{x}h\left(x\right)y-{\▿}_{x}g\left(x\right)(z+w)=0\\ L_y\≡h\left(x\right)=0\\ L_z\≡g\left(x\right)-l-\underset{\‾}{g}=0\\ L_w\≡g\left(x\right)+u-\stackrel{\‾}{g}=0\end{array}\right.;$

其中，扰动互补条件为 $\{\begin{array}{c}{L}_l^{\mathrm{\μ}}=LZe-\mathrm{\μ}e=0\\ L_u^{\mathrm{\μ}}=UWe+\mathrm{\μ}e=0\end{array},$

目标函数： $f(\·)=\underset{i\∈S_{GN}}{\Σ}{P}_{Gi}-\underset{i\∈S_{DN}}{\Σ}{P}_{Di}$

变量：x＝(P_G，Q_R，δ，V)，其中，h(x)为潮流方程(1)式；为潮流方程对状态变量的一阶偏导；为潮流方程对状态变量的二阶偏导；g(x)为不等式方程组；为不等式方程组对状态变量的一阶偏导；为不等式方程组对状态变量的二阶偏导；L表示以l_i为对角线的矩阵；Z表示以z_i为对角线的矩阵；W表示以w_i对角线的矩阵；是目标函数对状态变量的一阶导数；目标函数对状态变量的二阶导数；P_Gi：第i台机组有功出力；P_Di：第i节点有功负荷；S_GN：发电机集合；S_DN：负荷集合；c表示单位1的一维列向量；N为最大迭代次数，根据所运行的省地县电网大小所确定的有限循环次数；x＝(P_G，Q_R，δ，V)，P_G表示发电机节点的控制变量；Q_R表示PV节点的控制变量；δ表示所有节点的电压相角的状态变量；V表示所有节点的电压幅值状态变量；

(1-2-3)计算得出Δx、Δy、Δl、Δu、Δz、Δw，利用 $\left\{\begin{array}{c}{\mathrm{\&alpha;}}_p=0.9995\min \{\underset{i}{\min }(\frac{-l_i}{{\mathrm{\&Delta;l}}_i},{\mathrm{\&Delta;l}}_i<0;\frac{-u_i}{{\mathrm{\&Delta;u}}_i},{\mathrm{\&Delta;u}}_i<0),1\}\\ {\mathrm{\&alpha;}}_d=0.9995\min \{\underset{i}{\min }(\frac{-z_i}{{\mathrm{\&Delta;z}}_i},{\mathrm{\&Delta;z}}_i<0;\frac{-w_i}{{\mathrm{\&Delta;w}}_i},{\mathrm{\&Delta;w}}_i>0),1\}\end{array}\right.,i=1,2,\mathrm{...},r,$ 计算迭代步长α_p和α_d；

(1-2-4)利用 $\left\{\begin{array}{c}{x}^{(k+1)}=x^{\left(k\right)}+{\mathrm{\&alpha;}}_p\mathrm{\&Delta;}x\\ l^{(k+1)}=l^{\left(k\right)}+{\mathrm{\&alpha;}}_p\mathrm{\&Delta;}l\\ u^{(k+1)}=u^{\left(k\right)}+{\mathrm{\&alpha;}}_p\mathrm{\&Delta;}u\end{array}\right.$ 和 $\left\{\begin{array}{c}{y}^{(k+1)}=y^{\left(k\right)}+{\mathrm{\&alpha;}}_d\mathrm{\&Delta;}y\\ z_l^{(k+1)}=z^{\left(k\right)}+{\mathrm{\&alpha;}}_d\mathrm{\&Delta;}z\\ w^{(k+1)}=w^{\left(k\right)}+{\mathrm{\&alpha;}}_d\mathrm{\&Delta;}w\end{array}\right.$ 计算x^(k+1)，l^(k+1)，u^(k+1)，y^{(k +1)}，z_l ^(k+1)和w^(k+1)；使k₁的值增加1，转入步骤(1-2-1)；

(1-3)若max|(ΔP_i，ΔQ_i)|＜ε，则将省地县数据文件进行Ward网络等值，得到网架数据和系统运行数据；

(1-4)重复步骤(1-1-1)至(1-1-4)计算网架数据和系统运行数据是否收敛，若k₁＞N，则用户控制计算机选择其它时间尺度的网架数据和系统运行数据并重复步骤(1-1-1)至(1-3)进行省地县潮流计算；

(1-5)如果潮流计算max|(ΔP_i，ΔQ_i)|＜ε，则从生成的等值网络数据文件中解析出多时间尺度的省地县一体化无功优化所需的网架数据及系统运行数据；

(1-6)对解析出的多时间尺度的省地县一体化无功优化所需的网架数据及运行数据进行联合优化：

(1-6-1)读取网架数据和运行数据，设定迭代次数k₂＝0，最大迭代次数限制为k_max＝M2，设定节点电压V_i＝1和节点相角δ_i＝0，设定松弛变量l_i＝1，u_i＝1，拉格朗日乘子z_i＝1，w_i＝-0.5，y_i＝0；其中，中心参数σ∈(0，1)，根据网架参数生成导纳矩阵Y_B；

(1-6-2)当k₂＜M2时，计算互补间隙Gap＝l^Tz-u^Tw，如果互补间隙Gap＜ε，则最优无功计算收敛，生成包含电网各个节点无功功率分布、变电站无功设备投切数量以及变压器分接头位置的结果文件，其中收敛精度ε＝10^-6；若Gap≥ε，则利用μ＝σ(l^Tz-u^Tw)/(2r)计算扰动因子μ，并计算修正方程组

$\{\begin{array}{c}-\&lsqb;{\&dtri;}_x^{2}f\left(x\right)-{\&dtri;}_x^{2}h\left(x\right)y-{\&dtri;}_x^{2}g\left(x\right)(z+w)\&rsqb;\mathrm{\&Delta;}x+{\&dtri;}_{x}h\left(x\right)\mathrm{\&Delta;}y+{\&dtri;}_{x}g\left(x\right)(\mathrm{\&Delta;}z+\mathrm{\&Delta;}w)=L_x\\ {\&dtri;}_x{h}^T\left(x\right)\mathrm{\&Delta;}x=-L_y\\ {\&dtri;}_x{g}^T\left(x\right)\mathrm{\&Delta;}x-\mathrm{\&Delta;}x=-L_z\\ {\&dtri;}_x{g}^T\left(x\right)\mathrm{\&Delta;}x+\mathrm{\&Delta;}u=-L_w\\ Z\mathrm{\&Delta;}l+L\mathrm{\&Delta;}z=-L_l^{\mathrm{\&mu;}}\\ W\mathrm{\&Delta;}u+U\mathrm{\&Delta;}w=-L_u^{\mathrm{\&mu;}}\end{array};$

$\left\{\begin{array}{c}{L}_x\&equiv;{\&dtri;}_{x}f\left(x\right)-{\&dtri;}_{x}h\left(x\right)y-{\&dtri;}_{x}g\left(x\right)(z+w)=0\\ L_y\&equiv;h\left(x\right)=0\\ L_z\&equiv;g\left(x\right)-l-\underset{\&OverBar;}{g}=0\\ L_w\&equiv;g\left(x\right)+u-\stackrel{\&OverBar;}{g}=0\end{array}\right.$

扰动互补条件 $\{\begin{array}{c}{L}_l^{\mathrm{\&mu;}}=LZe-\mathrm{\&mu;}e=0\\ L_u^{\mathrm{\&mu;}}=UWe+\mathrm{\&mu;}e=0\end{array},$

$\begin{array}{c}{P}_{ij}^{set}={\left(V_i^{set}\right)}^2{Y}_{ij}^{set}{\mathrm{cos\&alpha;}}_{ij}^{set}-V_i^{set}{Y}_j^{set}{\mathrm{cos\&delta;}}_{ij}^{set},\\ {\left(I_{ij}^{set}\right)}^2={\left(V_i^{set}\right)}^2{\left(Y_{ij}^{set}\right)}^2+{\left(V_j^{set}\right)}^2{\left(Y_{ij}^{set}\right)}^2-2V_i^{set}{\left(Y_{ij}^{set}\right)}^2{V}_j^{set}\cos ({\mathrm{\&delta;}}_i^{set}+{\mathrm{\&alpha;}}_{ij}^{set})\cos ({\mathrm{\&delta;}}_j^{set}+{\mathrm{\&alpha;}}_{ij}^{set})\\ -2V_i^{set}{\left(Y_{ij}^{set}\right)}^2{V}_j^{set}\sin ({\mathrm{\&delta;}}_i^{set}+{\mathrm{\&alpha;}}_{ij}^{set})\sin ({\mathrm{\&delta;}}_j^{set}+{\mathrm{\&alpha;}}_{ij}^{set})\end{array}$

其中，r为不等式约束的个数，

$h\left(x\right)=\left\{\begin{array}{cc}{P}_{Gi}^{set}-P_{Di}^{set}-\underset{j\&Element;S_N^{set}}{\&Sigma;}{V}_i^{set}{Y}_{ij}^{set}{V}_j^{set}{\mathrm{cos\&delta;}}_{ij}^{set}=0& i\&Element;S_N^{set}\\ Q_{Gi}^{set}-Q_{Di}^{set}-\underset{j\&Element;S_N^{set}}{\&Sigma;}{V}_i^{set}{Y}_{ij}^{set}{V}_j^{set}{\mathrm{sin\&delta;}}_{ij}^{set}=0& i\&Element;S_N^{set}\\ (k_{ij}^{set1}-k_{ij1}^{set1})(k_{ij}^{set1}-k_{ij2}^{set1})\mathrm{...}(k_{ij}^{set1}-k_{ijm}^{set1})=0& i,j\&Element;S_K^{set}\\ (Q_{Ci}^{set}-Q_{Ci1}^{set})(Q_{Ci}^{set}-Q_{Ci2}^{set})\mathrm{...}(Q_{Ci}^{set}-Q_{Cim}^{set})=0& i\&Element;S_T^{set}\end{array}\right.$

目标函数： $f(\&CenterDot;)=\left\{\begin{array}{c}\underset{i\&Element;S_{GN}^{set}}{\&Sigma;}{P}_{Gi}^{set}-\underset{i\&Element;S_{DN}^{set}}{\&Sigma;}{P}_{Di}^{set}\\ \underset{i\&Element;S_{GN}^{set}}{\&Sigma;}{Q}_{Gi}^{set}-\underset{i\&Element;S_{DN}^{set}}{\&Sigma;}{Q}_{Di}^{set}\end{array}\right.$

变量：x＝(P_G，Q_R，δ，V)

其中，h(x)为潮流方程式；为潮流方程对状态变量的一阶偏导；为潮流方程对状态变量的二阶偏导；g(x)为不等式方程组；为不等式方程组对状态变量的一阶偏导；为不等式方程组对状态变量的二阶偏导；L表示以l_i为对角线的矩阵；Z表示以z_i为对角线的矩阵；W表示以w_i对角线的矩阵；是目标函数对状态变量的一阶导数；目标函数对状态变量的二阶导数；α_ij：节点导纳矩阵对应元素的夹角；省地县第i台机组有功出力；省地县第i节点有功负荷；省地县第i台机组无功出力；省地县第i节点无功负荷；省地县发电机集合；省地县负荷集合；为省地该变压器抽头的各档位对应值，m为其总档位，为省地可调变压器组集合；为该可调电容组投切组数的对应无功出力值，m为其总组数，为省地县可调电容组集合；c表示单位1的一维列向量；M2为最大迭代次数，根据所运行的省地县电网大小所确定的有限循环次数；上标set∈{s，x，d}的含义分别为：s为省，d为地，x为县；由于县级变压器不可调，因此调压器分接头调节只针对省地所属的变压器，上标set1∈{s，d}的含义分别为：s为省，d为地；x＝(P_G，Q_R，δ，V)，P_G表示发电机节点的控制变量；Q_R表示PV节点的控制变量；δ表示所有节点的电压相角的状态变量；V表示所有节点的电压幅值状态变量；

(1-6-3)计算得出Δx、Δy、Δl、Δu、Δz、Δw，利用 $\left\{\begin{array}{c}{\mathrm{\&alpha;}}_p=0.9995\min \{\underset{i}{\min }(\frac{-l_i}{{\mathrm{\&Delta;l}}_i},{\mathrm{\&Delta;l}}_i<0;\frac{-u_i}{{\mathrm{\&Delta;u}}_i},{\mathrm{\&Delta;u}}_i<0),1\}\\ {\mathrm{\&alpha;}}_d=0.9995\min \{\underset{i}{\min }(\frac{-z_i}{{\mathrm{\&Delta;z}}_i}{\mathrm{\&Delta;z}}_i<0;\frac{-w_i}{{\mathrm{\&Delta;w}}_i},{\mathrm{\&Delta;w}}_i>0),1\}\end{array}\right.,i=1,2,\mathrm{...},r,$ 计算迭代步长α_p和α_d；

(1-6-4)利用 $\left\{\begin{array}{c}{x}^{(k+1)}=x^{\left(k\right)}+{\mathrm{\&alpha;}}_p\mathrm{\&Delta;}x\\ l^{(k+1)}=l^{\left(k\right)}+{\mathrm{\&alpha;}}_p\mathrm{\&Delta;}l\\ u^{(k+1)}=u^{\left(k\right)}+{\mathrm{\&alpha;}}_p\mathrm{\&Delta;}u\end{array}\right.$ 和 $\left\{\begin{array}{c}{y}^{(k+1)}=y^{\left(k\right)}+{\mathrm{\&alpha;}}_d\mathrm{\&Delta;}y\\ z_l^{(k+1)}=z^{\left(k\right)}+{\mathrm{\&alpha;}}_d\mathrm{\&Delta;}z\\ w^{(k+1)}=w^{\left(k\right)}+{\mathrm{\&alpha;}}_d\mathrm{\&Delta;}w\end{array}\right.$ 计算x^(k+1)，l^(k+1)，u^(k+1)，y^{(k +1)}，z_l(^k+1)和w^(k+1)；使k₂值增加1，转入步骤(1-6-2)；

(1-7)若多时间尺度的省地县一体化无功优化计算互补间隙Gap＜ε，则无功优化收敛，根据无功优化计算对主变分接开关档位调节、无功补偿器投切量以及部分变电站增加电容器和电抗器详细情况更新初始的数据文件并生成新的省地县数据文件；

(1-8)如果迭代次数k₂＞M2意味着多时间尺度的省地县一体化无功优化计算不收敛或生成的省地县数据文件校验结果不收敛，则调整节点电压不等式约束的上下限，转入步骤(1-6-1)；

(1-9)根据无功优化模块计算收敛后生成的结果文件，结果文件包括整个系统的无功功率分布、各节点的无功补偿量、变压器变比分接头位置、无功补偿装置投切量、节点电压分布、系统的有功损耗、系统的无功损耗、省与地之间的无功功率交换量、地与县之间的无功功率交换量、发电厂的有功出力和无功出力，对数据进行分类形成标准的电网分析报表；

(1-10)计算机根据电网分析报表中的分析结果及运行方式调整策略报告，调节主变分接开关档位，调节电网中每个无功补偿器投切量并调节电网中部分变电站的电容器和\或电抗器的数量，确定每座变电站的变压器抽头档位、电容器投切组数、电网中发电机节点的本地发电机机端电压和功率。

本发明的多时间尺度的省地县一体化无功优化方法是一种考虑了年规律、季规律、月规律、日规律、超短时固有惯性规律和随机成份的省地县一体化无功优化方法。能够实现多时间尺度的省地县一体化无功优化，给出多时间尺度连续的无功电压优化配置策略，实现电力系统长、中、短期各电压等级间的无功电压优化运行，提高全网的电压合格率，保障电网的持续的安全稳定运行，为制定电网的运行方式、合理且优化的利用发电资源提供理论依据，能够有效的确保电网安全稳定持续运行。

作为优选，时间尺度为5个，分别为年度、季度、月度、天、滚动超短期；滚动超短期为15分钟。

作为优选，

所述网架数据包含省地县整个地区或部分区域的网络线路参数，接地支路参数和变压器参数，网络线路参数包括电抗值x_lineij、电阻值r_lineij和导纳值b_lineij，接地支路参数包括电导值G_i和电纳值B_i，变压器参数包括电抗值x_lineij、电阻值r_lineij和变比k_lineij；

系统运行参数模块包含省地县的整个网络或部分区域的运行参数和PV节点参数，运行参数包括发电机节点运行参数P_Gi、Q_Gi，和负荷节点的运行参数P_Di、Q_Di，PV节点参数包括节点电压V_i和无功功率上下限Q_i、

作为优选，M为3至20，N为3至20，M2为10至100。

作为优选，结果文件包括整个系统的无功功率分布、各节点的无功补偿量、变压器变比分接头位置、无功补偿装置投切量、节点电压分布、系统的有功损耗、系统的无功损耗、省与地之间的无功功率交换量、地与县之间的无功功率交换量、发电厂的有功出力和无功出力。

与现有技术方案相比，本发明的有益效果是：

针对目前电网无功优化管理分散、繁杂的特点，设计并实现了一种先进的多时间尺度的省地县一体化无功优化计算放法，计算过程可以全网联立计算，也可以通过边界条件并行计算，计算方式灵活。该方法将省地县高中低网络无功优化计算、分析集成于一体，利用科学的无功优化计算、比较、分析手段，为实现生产决策提供技术支持。方法主要分三个等级进行无功优化，主要为：单独无功优化、分级无功优化、全局无功优化。单独无功优化针对某地区的负荷水平以及电压质量进行跟踪补偿，实现就地补偿和就地平衡，降低线损，提高补偿效率。分级无功优化针对不同的网区电压等级以及用户负荷类型的，需要无功补偿的方式不同。全局无功优化是针对单独无功优化以及分级无功优化不合理方案可能导致局部无功供过于求或者整体无功不平衡等问题进行无功优化，再加上该方法还考虑多时段时间尺度，更能够充分性地考虑了电网复杂的电网特性规律进行无功优化，能够使整个地区无功设备最优利用，系统无功达到平衡，无功补偿设备费用减少，同时也减小了电网的网络损耗，保证各电压等级间的电压质量，满足国家有关法律法规相关技术标准。解决了传统无功优化的多方面不足，在国内外均没有先例。

附图说明

图1为本发明的实施例的一种流程图。

具体实施方式

下面参照附图，详细叙述本发明的具体实施方案。

如图1所示实施例是一种联合省地县电力系统多时间尺度的省地县一体化无功优化方法，包括如下步骤：

步骤100选择要优化的电网范围，对应其电网的网架数据和电网的运行数据，计算机中设有2011年年度、2011年夏季、2011年5月、2011年5月15日，用户控制计算机选择时间尺度为2011年年度的网架数据和运行数据进行省地县潮流计算：

步骤101，读取网架数据和运行数据，网架数据包括电网网架参数，运行数据包括系统运行的参数，设定V_i＝1和δ_i＝0，设定迭代次数k＝0，最大迭代次数k_max＝M，根据网架参数生成导纳矩阵Y_B；导纳矩阵Y_B中包含电网的网架结构关系；

步骤102，利用潮流方程，计算功率不平衡量ΔP_i和ΔQ_i：

步骤103，校验功率不平衡量ΔP_i和ΔQ_i的值，若max|(ΔP_i，ΔQ_i)|＜ε，则收敛，然后进行网络等值，其中收敛精度ε＝10^-6；

若max|(ΔP_i，ΔQ_i)|≥ε利用 $\left[\begin{array}{cc}\frac{\&part;\mathrm{\&Delta;}P}{\&part;\mathrm{\&delta;}}& \frac{\&part;\mathrm{\&Delta;}P}{\&part;V}\\ \frac{\&part;\mathrm{\&Delta;}Q}{\&part;\mathrm{\&delta;}}& \frac{\&part;\mathrm{\&Delta;}Q}{\&part;V}\end{array}\right]\left[\begin{array}{c}\mathrm{\&Delta;}\mathrm{\&delta;}\\ \mathrm{\&Delta;}V\end{array}\right]=-\left[\begin{array}{c}\mathrm{\&Delta;}P\\ \mathrm{\&Delta;}Q\end{array}\right]$ 计算不平衡量ΔV和Δδ；

步骤104，利用 $\left\{\begin{array}{c}\mathrm{\&Delta;}V=V^{(k+1)}-V^{\left(k\right)}\\ \mathrm{\&Delta;}\mathrm{\&delta;}={\mathrm{\&delta;}}^{(k+1)}-{\mathrm{\&delta;}}^{\left(k\right)}\end{array}\right.$ 求解V^(k+1)和δ^(k+1)，使k值增加1，当max|(ΔP_i，ΔQ_i)|≥ε时，重复步骤102；

步骤200，若k＞20，则通过最优潮流对网架参数、发电机的有功功率PG_i、无功功率QG_i出力大小进行调整计算：

步骤201，重新读取网架数据和运行数据，设定迭代次数k₁＝0，最大迭代次数限制为k_1max＝20，设置节点电压V_i＝1和节点相角δ_i＝0，设置松弛变量l_i＝1，u_i＝1，拉格朗日乘子z_i＝1，w_i＝-0.5，y_i＝0；σ∈(0，1)称为中心参数，根据网架参数生成导纳矩阵Y_B；

步骤202，当k₁＜20时，计算互补间隙Gap＝l^Tz-u^Tw，如果互补间隙Gap＜ε，则最优潮流收敛，然后进行网络等值，其中收敛精度ε＝10^-6；

若Gap≥ε，则利用μ＝σ(l^Tz-u^Tw)/(2r)计算扰动因子μ，并计算修正方程组 $\{\begin{array}{c}-\&lsqb;{\&dtri;}_x^{2}f\left(x\right)-{\&dtri;}_x^{2}h\left(x\right)y-{\&dtri;}_x^{2}g\left(x\right)(z+w)\&rsqb;\mathrm{\&Delta;}x+{\&dtri;}_{x}h\left(x\right)\mathrm{\&Delta;}y+{\&dtri;}_{x}g\left(x\right)(\mathrm{\&Delta;}z+\mathrm{\&Delta;}w)=L_x\\ {\&dtri;}_x{h}^T\left(x\right)\mathrm{\&Delta;}x=-L_y\\ {\&dtri;}_x{g}^T\left(x\right)\mathrm{\&Delta;}x-\mathrm{\&Delta;}x=-L_z\\ {\&dtri;}_x{g}^T\left(x\right)\mathrm{\&Delta;}x+\mathrm{\&Delta;}u=-L_w\\ Z\mathrm{\&Delta;}l+L\mathrm{\&Delta;}z=-L_l^{\mathrm{\&mu;}}\\ W\mathrm{\&Delta;}u+U\mathrm{\&Delta;}w=-L_u^{\mathrm{\&mu;}}\end{array},$

$\left\{\begin{array}{c}{L}_x\&equiv;{\&dtri;}_{x}f\left(x\right)-{\&dtri;}_{x}h\left(x\right)y-{\&dtri;}_{x}g\left(x\right)(z+w)=0\\ L_y\&equiv;h\left(x\right)=0\\ L_z\&equiv;g\left(x\right)-l-\underset{\&OverBar;}{g}=0\\ L_w\&equiv;g\left(x\right)+u-\stackrel{\&OverBar;}{g}=0\end{array}\right.;$

步骤203，计算得出Δx、Δy、Δl、Δu、Δz、Δw，利用方程 $\{\begin{array}{c}{\mathrm{\&alpha;}}_p=0.9995\min \{\underset{i}{\min }(\frac{-l_i}{{\mathrm{\&Delta;l}}_i},{\mathrm{\&Delta;l}}_i<0;\frac{-u_i}{{\mathrm{\&Delta;u}}_i},{\mathrm{\&Delta;u}}_i<0),1\}\\ {\mathrm{\&alpha;}}_d=0.9995\min \{\underset{i}{\min }(\frac{-z_i}{{\mathrm{\&Delta;z}}_i}{\mathrm{\&Delta;z}}_i<0;\frac{-w_i}{{\mathrm{\&Delta;w}}_i},{\mathrm{\&Delta;w}}_i>0),1\}\end{array},i=1,2,\mathrm{...},r,$ 计算迭代步长α_p和α_d；

步骤204，利用 $\left\{\begin{array}{c}{x}^{(k+1)}=x^{\left(k\right)}+{\mathrm{\&alpha;}}_p\mathrm{\&Delta;}x\\ l^{(k+1)}=l^{\left(k\right)}+{\mathrm{\&alpha;}}_p\mathrm{\&Delta;}l\\ u^{(k+1)}=u^{\left(k\right)}+{\mathrm{\&alpha;}}_p\mathrm{\&Delta;}u\end{array}\right.$ 和 $\left\{\begin{array}{c}{y}^{(k+1)}=y^{\left(k\right)}+{\mathrm{\&alpha;}}_d\mathrm{\&Delta;}y\\ z_l^{(k+1)}=z^{\left(k\right)}+{\mathrm{\&alpha;}}_d\mathrm{\&Delta;}z\\ w^{(k+1)}=w^{\left(k\right)}+{\mathrm{\&alpha;}}_d\mathrm{\&Delta;}w\end{array}\right.$ 计算x^(k+1)，l^(k+1)，u^(k+1)，y^(k+1)，z_l ^(k+1)和w^(k+1)；使k₁的值增加1，转入步骤201；

步骤300，若max|(ΔP_i，ΔQ_i)|＜ε，则将省地县数据文件进行Ward网络等值，得到网架数据和系统运行数据；

步骤400，重复步骤101至104计算网架数据和系统运行数据是否收敛，若k₁＞N，则用户控制计算机选择其它时间尺度的网架数据和系统运行数据并重复步骤101至300进行省地县潮流计算；

步骤500，如果潮流计算max|(ΔP_i，ΔQ_i)|＜ε，则从生成的等值网络数据文件中解析出多时间尺度的省地县一体化无功优化所需的网架数据及系统运行数据；

步骤600，对解析出的多时间尺度的省地县一体化无功优化所需的网架数据及运行数据进行联合优化：

步骤601，根据所得的网架结构参数和系统运行参数与数据库给定的变压器、电容器、电抗器参数相结合，生成多时间尺度的省地县一体化无功优化计算模块所需要的数据文件；

步骤602，读取网架数据和运行数据，设定迭代次数k₂＝0，最大迭代次数限制为k_max＝50，设定节点电压V_i＝1和节点相角δ_i＝0，设定松弛变量l_i＝1，u_i＝1，拉格朗日乘子z_i＝1，w_i＝-0.5，y_i＝0；其中，中心参数σ∈(0，1)，根据网架参数生成导纳矩阵Y_B；

步骤603，当k₂＜50时，计算互补间隙Gap＝l^Tz-u^Tw，如果互补间隙Gap＜ε，则最优无功计算收敛，生成包含电网各个节点无功功率分布、变电站无功设备投切数量以及变压器分接头位置的结果文件，其中收敛精度ε＝10^-6；若Gap≥ε，则利用μ＝σ(l^Tz-u^Tw)/(2r)计算扰动因子μ，并计算修正方程组

$\{\begin{array}{c}-\&lsqb;{\&dtri;}_x^{2}f\left(x\right)-{\&dtri;}_x^{2}h\left(x\right)y-{\&dtri;}_x^{2}g\left(x\right)(z+w)\&rsqb;\mathrm{\&Delta;}x+{\&dtri;}_{x}h\left(x\right)\mathrm{\&Delta;}y+{\&dtri;}_{x}g\left(x\right)(\mathrm{\&Delta;}z+\mathrm{\&Delta;}w)=L_x\\ {\&dtri;}_x{h}^T\left(x\right)\mathrm{\&Delta;}x=-L_y\\ {\&dtri;}_x{g}^T\left(x\right)\mathrm{\&Delta;}x-\mathrm{\&Delta;}x=-L_z\\ {\&dtri;}_x{g}^T\left(x\right)\mathrm{\&Delta;}x+\mathrm{\&Delta;}u=-L_w\\ Z\mathrm{\&Delta;}l+L\mathrm{\&Delta;}z=-L_l^{\mathrm{\&mu;}}\\ W\mathrm{\&Delta;}u+U\mathrm{\&Delta;}w=-L_u^{\mathrm{\&mu;}}\end{array};$

$\left\{\begin{array}{c}{L}_x\&equiv;{\&dtri;}_{x}f\left(x\right)-{\&dtri;}_{x}h\left(x\right)y-{\&dtri;}_{x}g\left(x\right)(z+w)=0\\ L_y\&equiv;h\left(x\right)=0\\ L_z\&equiv;g\left(x\right)-l-\underset{\&OverBar;}{g}=0\\ L_w\&equiv;g\left(x\right)+u-\stackrel{\&OverBar;}{g}=0\end{array}\right.;$

步骤604，计算得出Δx、Δy、Δl、Δu、Δz、Δw，利用方程 $\left\{\begin{array}{c}{\mathrm{\&alpha;}}_p=0.9995\min \{\underset{i}{\min }(\frac{-l_i}{{\mathrm{\&Delta;l}}_i},{\mathrm{\&Delta;l}}_i<0;\frac{-u_i}{{\mathrm{\&Delta;u}}_i},{\mathrm{\&Delta;u}}_i<0),1\}\\ {\mathrm{\&alpha;}}_d=0.9995\min \{\underset{i}{\min }(\frac{-z_i}{{\mathrm{\&Delta;z}}_i}{\mathrm{\&Delta;z}}_i<0;\frac{-w_i}{{\mathrm{\&Delta;w}}_i},{\mathrm{\&Delta;w}}_i>0),1\}\end{array}\right.,i=1,2,\mathrm{...},r,$ 计算迭代步长α_p和α_d；

步骤605，利用 $\left\{\begin{array}{c}{x}^{(k+1)}=x^{\left(k\right)}+{\mathrm{\&alpha;}}_p\mathrm{\&Delta;}x\\ l^{(k+1)}=l^{\left(k\right)}+{\mathrm{\&alpha;}}_p\mathrm{\&Delta;}l\\ u^{(k+1)}=u^{\left(k\right)}+{\mathrm{\&alpha;}}_p\mathrm{\&Delta;}u\end{array}\right.$ 和 $\left\{\begin{array}{c}{y}^{(k+1)}=y^{\left(k\right)}+{\mathrm{\&alpha;}}_d\mathrm{\&Delta;}y\\ z_l^{(k+1)}=z^{\left(k\right)}+{\mathrm{\&alpha;}}_d\mathrm{\&Delta;}z\\ w^{(k+1)}=w^{\left(k\right)}+{\mathrm{\&alpha;}}_d\mathrm{\&Delta;}w\end{array}\right.$ 计算x^(k+1)，l^(k+1)，u^(k+1)，y^(k+1)，z_l ^(k+1)和w^(k+1)；使k₂值增加1，转入步骤603；

步骤700，若多时间尺度的省地县一体化无功优化计算互补间隙Gap＜ε，则无功优化收敛，根据无功优化计算对主变分接开关档位调节、无功补偿器投切量以及部分变电站增加电容器和电抗器详细情况更新初始的数据文件并生成新的省地县数据文件；

步骤800，如果迭代次数k₂＞50意味着多时间尺度的省地县一体化无功优化计算不收敛或生成的省地县数据文件校验结果不收敛，则调整节点电压不等式约束的上下限，转入步骤602；

步骤900，根据无功优化模块计算收敛后生成的结果文件，对数据进行分类形成标准的电网分析报表；

步骤1000，计算机根据电网分析报表中的分析结果及运行方式调整策略报告，调节主变分接开关档位，调节电网中每个无功补偿器投切量并调节电网中部分变电站的电容器和\或电抗器的数量，确定每座变电站的变压器抽头档位、电容器投切组数、电网中发电机节点的本地发电机机端电压和功率。

步骤1100，将多时间尺度的省地县一体化无功优化计算模块计算前和计算后的无功电压控制设备的状态值进行统计形成无功电压分析报告，包括：500kV变压器抽头位置及对应电压值、500kV变电站无功补偿容量汇总表、220kV变电站无功补偿容量汇总表、电网统调机组出力汇总表、优化前后网损统计情况；110kV变电站无功补偿容量汇总表、电网统调机组出力汇总表、优化前后网损统计情况；35kV变电站无功补偿容量汇总表、优化前后网损统计情况；10kV优化前后网损统计情况；

步骤1200，通过分析多时间尺度的省地县一体化无功优化计算模块的计算结果给出电网省地县无功一体化配置的优化策略，用于指导实际运行方式的调整，包括：变压器各电压等级抽头档位调整、当前无功补偿配置情况下的电容器和电抗器投切量、部分变电站增加电容器和电抗器配置的建议，经过这些策略的调整，会使电网的有功功率损耗以及无功功率损耗降低，电网的电压质量提高，系统的稳定性加强。

本发明的多时间尺度的省地县一体化无功优化计算是在给定系统网络结构和系统负荷的情况下，保证满足用户用电需求的前提下，优化电力系统的静态运行条件，通过调节控制变量使指定的目标函数值达到最小，同时满足系统对控制变量、状态变量及变量函数的物理限制和运行限制。调节控制变量包括发电机的有功无功出力，变压器分接头，电容器组分组投切等。通过无功优化不仅可以保证全网的电压质量，而且能取得可观的经济效益，能实现电能质量、系统运行的安全性和经济性。

本发明中的控制变量参数设置包括：有功出力、无功出力、变压器抽头、电容电抗器投切组数。其中：有功出力设置为按照要求将一些电厂的有功出力设为可调；无功出力要求有功有限制的PV节点、电网外等值节点和电网的所有平衡节点的无功出力均可调；变压器抽头均可调；电网500kV、220kV、110kV、35kV、10kV变电站电容电抗器均可调。所有变电站的电容器最大允许补偿容量为现有配置补偿容量的80％，所有变电站的电抗器最大允许补偿容量均为现有电抗器配置容量。

其约束条件参数设置包括：无功负荷模型、节点电压。其中：无功负荷模型全部加入约束限制；节点电压限制包括500kV、220kV、110kV、35kV、10kV母线电压，发电机端电压，限制范围如表1所示。

表1节点电压限制范围：

表2多时间尺度的省地县一体化无功优化前后某电网500kV主变分接头位置：

本发明实施例的优化前后某电网500kV主变分接头位置如表2所示。通过联合优化，找出满足电网网损最小的各变压器分接头。

表3电网某多时间尺度的省地县一体化无功优化前后大方式网损率：

表4电网某多时间尺度的省地县一体化无功优化前后小方式网损率：

根据多时间尺度省地县运行方式无功优化联调大方式下优化前后的供电量、网损电量，通过网损率计算公式可计算出电网的网损率。某电网多时间尺度的省地县一体化无功优化前后大方式网损率情况如表3所示，优化后有功网损比优化前少5.89MW，网损率降低0.0353％。某电网多时间尺度的省地县一体化无功优化前后小方式网损率情况如表4所示，优化后有功网损比优化前少4.26MW，网损率降低0.0772％。

应理解，本实施例仅用于说明本发明而不用于限制本发明的范围。此外应理解，在阅读了本发明讲授的内容之后，本领域技术人员可以对本发明作各种改动或修改，这些等价形式同样落于本申请所附权利要求书所限定的范围。

Claims (5)

1.一种电力系统多时间尺度的省地县一体化无功优化方法，其特征是，包括下列步骤：

(1-1)选择要优化的电网范围，获得对应的电网的网架数据和电网的运行数据，计算机中设有若干个时间尺度，用户控制计算机选择至少一个1个时间尺度的网架数据和运行数据进行省地县潮流计算：

(1-1-1)读取网架数据和运行数据，网架数据包括网架参数，运行数据包括系统运行的参数，设定V_i＝1和δ_i＝0，V_i为第i个节点的电压幅值，δ_i为第i个节点的电压相角，i为系统节点编号，设定迭代次数k＝0，最大迭代次数k_max＝M，根据网架参数生成导纳矩阵Y_B；导纳矩阵Y_B中包含电网的网架结构关系；

(1-1-2)利用潮流方程，计算功率不平衡量ΔP_i和ΔQ_i；

(1-1-3)校验功率不平衡量ΔP_i和ΔQ_i的值，若max|(ΔP_i，ΔQ_i)|＜ε，则收敛，然后进行网络等值，其中收敛精度ε＝10^-6；

若max|(ΔP_i，ΔQ_i)|≥ε，利用计算不平衡量ΔV和Δδ；

其中，为有功功率偏量对节点电压相角的偏导；为有功功率偏量对电压的偏导；为无功功率偏量对节点电压相角的偏导；为无功功率偏量对电压的偏导；Δδ为节点电压相角增量；ΔV为节点电压幅值增量，ΔP为一维列向量：ΔP＝(ΔP₁，ΔP₂，…ΔP_i)，i∈S_N，ΔP_i为第i节点的有功功率偏量；ΔQ为一维列向量：ΔQ＝(ΔQ₁，ΔQ₂，…ΔQ_i)，i∈S_N，ΔQ_i为第i节点的无功功率偏量；S_N为电网节点集合；

(1-1-4)利用求解V^(k+1)和δ^(k+1)，V^(k+1)为第k+1次迭代的电压幅值，V^(k)为第k次迭代的电压幅值，δ^(k+1)为第k+1次迭代的电压相角，δ^(k)为第k次迭代的电压相角，使k值增加1，当max|(ΔP_i，ΔQ_i)|≥ε时，重复步骤(1-1-2)；

(1-2)若k＞M，则通过最优潮流对网架参数、发电机的有功功率P_Gi、无功功率Q_Gi出力大小进行调整计算：

(1-2-1)重新读取网架数据和运行数据，设定迭代次数k₁＝0，最大迭代次数限制为k_1max＝N，设置节点电压V_i＝1和节点相角δ_i＝0，设置松弛变量l_n＝1，u_n＝1，拉格朗日乘子z_n＝1，w_n＝-0.5，y_n＝0；σ∈(0，1)称为中心参数，根据网架参数生成导纳矩阵Y_B；

(1-2-2)当k₁＜N时，计算互补间隙Gap＝l^Tz-u^Tw，如果互补间隙Gap＜ε，则最优潮流收敛，然后进行网络等值，其中收敛精度ε＝10^-6，r为常数，S_r＝1，2，…，r；

其中，l为一维列向量，l^T是l的转置，l_n∈(l₁，l₂，…，l_r)，n∈S_r，l_n为第n个下界松弛变量；u为一维列向量，u^T是u的转置，u_n∈(u₁，u₂，…u_r)，n∈S_r，u_n为第n个上界松弛变量；z为一维列向量：z_n∈(z₁，z₂，…z_r)，n∈S_r，z_n为第n个下界不等式约束的拉格朗日乘子；w为一维列向量：w_n∈(w₁，w₂，…w_r)，n∈S_r，w_n为第n个上界不等式约束的拉格朗日乘子；

若Gap≥ε，则利用μ＝σ(l^Tz-u^Tw)/(2r)计算扰动因子μ，并计算修正方程组

$\left\{\begin{array}{c}{L}_x\&equiv;{\&dtri;}_{x}f\left(x\right)-{\&dtri;}_{x}h\left(x\right)y-{\&dtri;}_{x}g\left(x\right)(z+w)=0\\ L_y\&equiv;h\left(x\right)=0\\ L_z\&equiv;g\left(x\right)=-\underset{\&OverBar;}{g}=0\\ L_w\&equiv;g\left(x\right)+u-\stackrel{\&OverBar;}{g}=0\end{array}\right.;$

变量：x＝(P_G，Q_R，δ，V)；Δx，Δy，Δl，Δu，Δz，Δw分别表示变量x，y，l，u，z，w的修正量；y为等式约束h(x)对应的拉格朗日乘子变量；h(x)为潮流方程式；为潮流方程对状态变量的一阶偏导；为潮流方程对状态变量的一阶偏导的转置；为潮流方程对状态变量的二阶偏导；g(x)为不等式方程组；为不等式方程组对状态变量的一阶偏导；为不等式方程组对状态变量的一阶偏导的转置；为不等式方程组对状态变量的二阶偏导；是目标函数对状态变量的一阶导数；目标函数对状态变量的二阶导数；L表示以l为对角线的矩阵；Z表示以z为对角线的矩阵；W表示以w对角线的矩阵；U表示以u为对角线的矩阵；L_x，L_y，L_z，L_w，和分别表示对变量x，y，z，w，l，u求偏导后在给定点处的数值；g和分别表示不等式约束的下界和上界；

(1-2-3)计算得出Δx、Δy、Δl、Δu、Δz、Δw，利用方程计算迭代步长α_p和α_d，Δl_n、Δu_n、Δz_n、Δw_n分别表示变量l_n、u_n、z_n、w_n的增量；

(1-2-4)利用和计算 和和分别表示x，l，u，y，z，w第k₁+1次的迭代值，k₁表示第k₁次迭代；使k₁的值增加1，转入步骤(1-2-1)；

(1-3)若max|(ΔP_i，ΔQ_i)|＜ε，则将省地县数据文件进行Ward网络等值，得到网架数据和系统运行数据；

(1-4)重复步骤(1-1-1)至(1-1-4)计算网架数据和系统运行数据是否收敛，若k₁＞N，则用户控制计算机选择其它时间尺度的网架数据和系统运行数据并重复步骤(1-1-1)至(1-3)进行省地县潮流计算；

(1-5)如果潮流计算max|(ΔP_i，ΔQ_i)|＜ε，则从生成的等值网络数据文件中解析出多时间尺度的省地县一体化无功优化所需的网架数据及系统运行数据；

(1-6)对解析出的多时间尺度的省地县一体化无功优化所需的网架数据及运行数据进行联合优化：

(1-6-1)读取网架数据和运行数据，设定迭代次数k₂＝0，最大迭代次数限制为k_2max＝M2，设定节点电压和节点相角 为第i₁个节点的电压幅值，为第i₁个节点的电压相角，i₁为省地县一体化无功优化系统节点编号，设置松弛变量l_n＝1，u_n＝1，拉格朗日乘子z_n＝1，w_n＝-0.5，y_n＝0；σ∈(0，1)称为中心参数，根据网架参数生成导纳矩阵Y_B；

(1-6-2)当k₂＜M2时，计算互补间隙Gap＝l^Tz-u^Tw，如果互补间隙Gap＜ε，则最优无功计算收敛，生成包含电网各个节点无功功率分布、变电站无功设备投切数量以及变压器分接头位置的结果文件，其中收敛精度ε＝10^-6；

若Gap≥ε，则利用μ＝σ(l^Tz-u^Tw)/(2r)计算扰动因子μ，并计算修正方程组

$\left\{\begin{array}{c}{L}_x\&equiv;{\&dtri;}_{x}f\left(x\right)-{\&dtri;}_{x}h\left(x\right)y-{\&dtri;}_{x}g\left(x\right)(z+w)=0\\ L_y\&equiv;h\left(x\right)=0\\ L_z\&equiv;g\left(x\right)=-\underset{\&OverBar;}{g}=0\\ L_w\&equiv;g\left(x\right)+u-\stackrel{\&OverBar;}{g}=0\end{array}\right.;$

扰动互补条件

$P_{i_1{j}_1}^{set}={\left(V_{i_1}^{set}\right)}^2{Y}_{i_1{j}_1}^{set}{\mathrm{cos\&alpha;}}_{i_1{j}_1}^{set}-V_{i_1}^{set}{Y}_{i_1{j}_1}^{set}{V}_{j_1}^{set}{\mathrm{cos\&delta;}}_{i_1{j}_1}^{set},$

$\begin{array}{c}{\left(I_{i_1{j}_1}^{set}\right)}^2={\left(V_{i_1}^{set}\right)}^2{\left(Y_{i_1{j}_1}^{set}\right)}^2+{\left(V_{j_1}^{set}\right)}^2{\left(Y_{i_1{j}_1}^{set}\right)}^2-2V_{i_1}^{set}{\left(Y_{i_1{j}_1}^{set}\right)}^2{V}_{j1}^{set}\cos ({\mathrm{\&delta;}}_{i_1}^{set}+{\mathrm{\&alpha;}}_{i_1{j}_1}^{set})\cos ({\mathrm{\&delta;}}_{j_1}^{set}+{\mathrm{\&alpha;}}_{i_1{j}_1}^{set})\\ -2V_{i_1}^{set}{\left(Y_{i_1{j}_1}^{set}\right)}^2{V}_{j_1}^{set}\sin ({\mathrm{\&delta;}}_{i_1}^{set}+{\mathrm{\&alpha;}}_{i_1{j}_1}^{set})\sin ({\mathrm{\&delta;}}_{j_1}^{set}+{\mathrm{\&alpha;}}_{i_1{j}_1}^{set})\end{array}$

$h\left(x\right)=\left\{\begin{array}{cc}{P}_{Gm}^{set}-P_{Dm}^{set}-{\mathrm{\&Sigma;V}}_{i_1}^{set}{Y}_{i_1{j}_1}^{set}{V}_{j_1}^{set}{\mathrm{cos\&delta;}}_{i_1{j}_1}^{set}=0& m\&Element;S_{GN}^{set};w\&Element;S_D^{set};i_1,j_1\&Element;S_N^{set}\\ Q_{Gm}^{set}-Q_{Dw}^{set}-{\mathrm{\&Sigma;V}}_{i_1}^{set}{Y}_{i_1{j}_1}^{set}{V}_{j_1}^{set}{\mathrm{sin\&delta;}}_{i_1{j}_1}^{set}=0& m\&Element;S_{GN}^{set};w\&Element;S_D^{set};i_1,j_1\&Element;S_N^{set}\\ (k_{uv}^{set1}-k_{uv1}^{set1})(k_{uv}^{set1}-k_{uv2}^{set1})...(k_{uv}^{set1}-k_{uvq}^{set1})=0& u,v\&Element;S_K^{set1}\\ (Q_{Ct}^{set}-Q_{Ct1}^{set})(Q_{Ct}^{set}-Q_{Ct2}^{set})...(Q_{Ct}^{set}-Q_{Ctp}^{set})=0& t\&Element;S_C^{set}\end{array}\right.$

目标函数：

变量：x＝(P_G，Q_R，δ，V)

其中，h(x)为潮流方程式；为潮流方程对状态变量的一阶偏导；为潮流方程对状态变量的一阶偏导的转置；为潮流方程对状态变量的二阶偏导；g(x)为不等式方程组；为不等式方程组对状态变量的一阶偏导；为不等式方程组对状态变量的一阶偏导的转置；为不等式方程组对状态变量的二阶偏导；是目标函数对状态变量的一阶导数；目标函数对状态变量的二阶导数；L表示以l为对角线的矩阵；Z表示以z为对角线的矩阵；W表示以w对角线的矩阵；U表示以u为对角线的矩阵；省地县第i₁个节点的电压幅值；省地县第j₁个节点的电压幅值；省地县电网节点集合；省地县第m台发电机有功功率；省地县第w节点有功负荷；省地县第m台机组无功出力；省地县第w节点无功负荷；省地县发电机集合；省地县负荷集合；为省地变压器分接头的各档位对应值，q为其总档位，为省地可调变压器组集合；为可调电容组投切组数的对应无功出力值，p为其总组数，为省地县可调电容组集合；为省地县支路有功集合；为省地县支路电流集合；为省地县线路i₁j₁的节点导纳幅值；为省地县线路i₁j₁的节点导纳角度；为省地县节点i₁的电压角度；为省地县节点j₁的电压角度；为省地县节点i₁和节点j₁电压角度之差；e表示单位1的一维列向量；M2为最大迭代次数，根据所运行的省地县电网大小所确定的有限循环次数；上标set∈{s，x，d}的含义分别为：s为省，d为地，x为县；由于县级变压器不可调，因此调压器分接头调节只针对省地所属的变压器，上标set1∈{s，d}的含义分别为：s为省，d为地；x＝(P_G，Q_R，δ，V)，P_G表示发电机节点的控制变量；Q_R表示PV节点的控制变量；δ表示所有节点的电压相角的状态变量；V表示所有节点的电压幅值状态变量；

(1-6-3)计算得出Δx、Δy、Δl、Δu、Δz、Δw，利用方程计算迭代步长α_p和α_d；

(1-6-4)利用和计算 和分别表示x，l，u，y，z，w第k₂+1次的迭代值，上标k₂表示第k₂次迭代；使k₂值增加1，转入步骤(1-6-2)；

(1-7)若多时间尺度的省地县一体化无功优化计算互补间隙Gap＜ε，则无功优化收敛，根据无功优化计算对主变分接开关档位调节、无功补偿器投切量以及部分变电站增加电容器和电抗器详细情况更新初始的数据文件并生成新的省地县数据文件；

(1-8)如果迭代次数k₂＞M2意味着多时间尺度的省地县一体化无功优化计算不收敛或生成的省地县数据文件校验结果不收敛，则调整节点电压不等式约束的上下限，转入步骤(1-6-1)；

(1-9)根据无功优化模块计算收敛后生成的结果文件，对数据进行分类形成标准的电网分析报表；

(1-10)计算机根据电网分析报表中的分析结果及运行方式调整策略报告，调节主变分接开关档位，调节电网中每个无功补偿器投切量并调节电网中部分变电站的电容器和\或电抗器的数量，确定每座变电站的变压器分接头档位、电容器投切组数、电网中发电机节点的本地发电机机端电压和功率。

2.根据权利要求1所述的电力系统多时间尺度的省地县一体化无功优化方法，其特征是，时间尺度为5个，分别为年度、季度、月度、天、滚动超短期；滚动超短期为15分钟。

3.根据权利要求1所述的电力系统多时间尺度的省地县一体化无功优化方法，其特征是，所述网架数据包含省地县整个地区或部分区域的网络线路参数，接地支路参数和变压器参数，网络线路参数包括电抗值电阻值和导纳值接地支路参数包括电导值和电纳值变压器参数包括电抗值电阻值和变比

系统运行参数模块包含省地县的整个网络或部分区域的运行参数和PV节点参数，运行参数包括发电机节点运行参数P_Gm、Q_Gm，和负荷节点的运行参数P_Dw、Q_Dw，PV节点参数包括节点电压和无功功率上下限

4.根据权利要求1所述的电力系统多时间尺度的省地县一体化无功优化方法，其特征是，M为3至20，N为3至20，M2为10至100。

5.根据权利要求1或2或3或4所述的电力系统多时间尺度的省地县一体化无功优化方法，其特征是，结果文件包括整个系统的无功功率分布、各节点的无功补偿量、变压器分接头位置、无功补偿装置投切量、节点电压分布、系统的有功损耗、系统的无功损耗、省与地之间的无功功率交换量、地与县之间的无功功率交换量、发电厂的有功功率和无功功率。