基于内点法最优潮流计算 PPT

基于功率—电流混合潮流约束的内点法最优潮流江全元,黄志光(浙江大学电气工程学院,浙江省杭州市310027)摘要:提出了一种基于直角坐标下功率—电流混合型潮流约束的最优潮流模型。

该模型对系统中的非零注入功率节点采用功率失配型潮流约束,而对零注入功率节点采用电流失配型潮流约束。

这种混合模型结合了功率和电流型潮流方程的优点:对于零注入功率节点,该模型具有电流型潮流方程一阶导数为常数、二阶导数为0的特点,从而使雅可比矩阵和海森矩阵更容易计算;对于非零注入功率节点,该模型也比电流型潮流方程更好处理。

该模型特别适合非线性预测—校正内点法的最优潮流,多个大规模算例测试证明该模型收敛性更好,计算效率更高,尤其适合求解含较多零注入功率节点的大规模电力系统最优潮流问题。

关键词:最优潮流;内点法;电流失配;功率失配;混合模型中图分类号:TM744收稿日期:2008211219;修回日期:2009203202。

教育部科学技术研究重点项目(107063);浙江省自然科学基金资助项目(R1080089)。

0　引言自从20世纪60年代Carpentier 提出最优潮流(O PF )问题以来,该问题受到越来越多的关注[127]。

Dommel 和Tinney 建立了OPF 问题的非线性规划模型。

实际应用中,OPF 有不同的表述形式,学者们也提出了不同的算法求解该问题。

近年来,内点法及其改进算法由于其突出的计算性能在OPF 问题中得到了广泛的应用[8214]。

O PF 问题中,网络潮流约束可表示为直角坐标或极坐标形式。

文献[15]比较了直角坐标和极坐标2种形式的OPF 模型,指出二者有相似的收敛特性和计算效率。

但在O PF 的内点算法中,直角坐标系的表达更简洁[16]。

在潮流计算时,系统潮流方程通常表示为功率失配形式,而电流失配型潮流方程突出了电力网络本身是线性、而节点注入是非线性的特点[17],一般来说更加适合求解负荷潮流问题[18]。

一种基于Karmarkar内点法的最优潮流算法
郝玉国;刘广一;于尔铿
【期刊名称】《中国电机工程学报》
【年(卷),期】1996(16)6
【摘要】以原－对偶内点算法（Ｋａｒｍａｒｋａｒ内点法的一种变形）为基本算法解算最优潮流问题，综合考虑非线性目标函数和约束条件，结合牛顿法最优潮流先进的稀疏矩阵技术，并且提出了一种新的原－对偶内点算法迭代步长选取原则和障碍参数修正策略。

算例表明本算法有较好的数值稳定性，优化结果精确，对不等式约束有较强的处理能力，显示了内点算法应用于大规模电力系统优化问题的良好前景。

【总页数】4页(P409-412)
【关键词】内点算法;最优潮流;电力系统
【作者】郝玉国;刘广一;于尔铿
【作者单位】电力科学研究院电网所
【正文语种】中文
【中图分类】TM744
【相关文献】
1.基于改进多中心校正解耦内点法的动态最优潮流并行算法 [J], 简金宝;杨林峰;全然
2.基于非线性多中心校正内点法的最优潮流算法 [J], 蔡广林;张勇军;任震
3.基于自动微分技术的内点法最优潮流算法 [J], 耿光超;江全元
4.基于退火粒子群和内点法的改进最优潮流算法 [J], 陈丽光;文波;聂一雄
5.基于预测校正内点法的最优潮流算法 [J], 陆子强
因版权原因，仅展示原文概要，查看原文内容请购买。

内点法最优潮流MATLAB算法clear;%clc;%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%初始化%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%数据加载n=input('请输入要计算的节点系统(5):')load Node5.txt;%节点数据load Branch5.txt;%支路数据load Generator5.txt;%发电机数据Node=Node5;Branch=Branch5;Generator=Generator5;%节点数据处理N=Node(:,1);%节点号Type=Node(:,2);%节点类型Uamp=Node(:,3);%节点电压幅值Dlta=Node(:,4);%节点电压相角Pd=Node(:,5);%节点负荷有功Qd=Node(:,6);%节点负荷无功Pg=Node(:,7);%节点出力有功Qg=Node(:,8);%节点出力无功Umax=Node(:,9);%节点电压幅值上限 Umin=Node(:,10);%节点电压幅值下限Bc=Node(:,11);%节点补偿电容电纳值 %支路数据处理Nbr=Branch(:,1);%支路号Nl=Branch(:,2);%支路首节点Nr=Branch(:,3);%支路末节点R=Branch(:,4);%支路电阻X=Branch(:,5);%支路电抗Z=R+1i*X;%支路阻抗=支路电阻+支路电抗 Bn=Branch(:,6);%支路对地电纳K=Branch(:,7);%支路变压器变比，0表示无变压器 Ptmax=Branch(:,8);%线路传输功率上限%发电机数据处理Ng=Generator(:,1);%发电机序号Nbus=Generator(:,2);%所在母线号Pumax=Generator(:,3);%发电机有功出力上界 Qumax=Generator(:,4);%发电机无功出力上界 Pumin=Generator(:,5);%发电机有功出力下界Qumin=Generator(:,6);%发电机无功出力下界a2=Generator(:,7);%燃料耗费曲线二次系数a1=Generator(:,8);%燃料耗费曲线一次系数a0=Generator(:,9);%燃料耗费曲线常数项%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%n=length(N);%节点个数ng=length(Ng);%发电机台数nbr=length(Nbr);%支路个数x=zeros(2*(ng+n),1);%控制变量+状态变量x(1:ng)=Pg(Nbus);x(ng+1:2*ng)=Qg(Nbus);x((2*ng+2):2:2*(ng+n))=Uamp; x((2*ng+1):2:2*(ng+n)-1)=Dlta; l=0.8*ones(2*ng+n+nbr,1);%松弛变量u=1.1*ones(2*ng+n+nbr,1);%松弛变量w=-1.5*ones(2*ng+n+nbr,1);%拉格朗日乘子z=ones(2*ng+n+nbr,1);%拉格朗日乘子y=zeros(2*n,1);%拉格朗日乘子y(1:2:2*n-1)=1e-3;y(2:2:2*n)=-1e-3;%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%计算不等式约束的上下限%%%%%%%%%%%%%%%%%%%%%%%%%gmingmin=zeros(2*ng+n+nbr,1);gmin(1:ng)=Pumin;gmin(ng+1:2*ng)=Qumin;gmin(2*ng+1:2*ng+n)=Umin;gmin(2*ng+n+1:2*ng+n+nbr)=-Ptmax; %gmaxgmax=zeros(2*ng+n+nbr,1);gmax(1:ng)=Pumax;gmax(ng+1:2*ng)=Qumax;gmax(2*ng+1:2*ng+n)=Umax;gmax(2*ng+n+1:2*ng+n+nbr)=Ptmax;%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%形成导纳矩阵%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%Y=zeros(n,n);%%%%%%%%%%%%%%%%%%%%计算非对角元素%%%%%%%%%%%%%%%%%%%%% for ii=1:nbr if K(ii)==0%非变压器支路Y(Nl(ii),Nr(ii))=-1/Z(ii);Y(Nr(ii),Nl(ii))=Y(Nl(ii),Nr(ii));else%变压器支路Y(Nl(ii),Nr(ii))=-1/Z(ii)/K(ii);Y(Nr(ii),Nl(ii))= Y(Nl(ii),Nr(ii));endend%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%计算对角元素%%%%%%%%%%%%%%%%%%%%%%for ii=1:n%将支路导纳加入到对角元素中for jj=1:nbrif K(jj)==0&&(Nl(jj)==ii||Nr(jj)==ii)%非变压器支路Y(ii,ii)=Y(ii,ii)+1/Z(jj);else if K(jj)~=0&&(Nl(jj)==ii||Nr(jj)==ii)%变压器支路Y(ii,ii)=Y(ii,ii)+1/Z(jj)/K(jj);endendendendfor ii=1:nbr%将对地电纳加入到对角元素中if K(ii)==0%非变压器支路Y(Nl(ii),Nl(ii))=Y(Nl(ii),Nl(ii))+1i*Bn(ii);Y(Nr(ii),Nr(ii))=Y(Nr(ii),Nr(ii))+1i*Bn(ii);else%变压器支路Y(Nr(ii),Nr(ii))=Y(Nr(ii),Nr(ii))+(K(ii)-1)/K(ii)/Z(ii);Y(Nl(ii),Nl(ii))=Y(Nl(ii),Nl(ii))+(1-K(ii))/K(ii)/K(ii)/Z(ii);endendfor ii=1:nY(ii,ii)=Y(ii,ii)+i*Bc(ii);endG=real(Y);%电导B=imag(Y);%电纳%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%k=0;%迭代次数Kmax=150;%最大迭代次数iteration=1e-4;%误差精度delta=0.08;Gap=(l'*z-u'*w)*ones(Kmax,1);%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%主程序%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% while k<50%计算互补间隙GapGap(k+1)=l'*z-u'*w;if Gap>iterationmiu=delta*Gap(k+1)/(2*(2*ng+n+nbr)); %%%%%%%%%%%%%%%%%%%%%%%%%%%%形成系数矩阵%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%相角差计算%%%%%%%%%%%%%%%%%%%%%%theta=zeros(n,n);for ii=1:nfor jj=1:ntheta(ii,jj)=Dlta(ii)-Dlta(jj);endend %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%1、等式约束雅克比矩阵%%%%%%%%%%%%%%%%pxh=zeros(2*(ng+n),2*n); %%%%%%%%%%%%%%%%%%%%%%%ah/aP%%%%%%%%%%%%%%%%%%%%%%%for ii=1:ngpxh(Ng(ii),2*Nbus(ii)-1)=1;end %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%ah/aQ%%%%%%%%%%%%%%%%%%%%%%%for ii=1:ngpxh(Ng(ii)+ng,2*Nbus(ii))=1;end %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%ah/ax%%%%%%%%%%%%%%%%%%%%%%%HH=zeros(n,n);JJ=zeros(n,n);NN=zeros(n,n);LL=zeros(n,n);for ii=1:nfor jj=1:nif ii~=jj%i!=j时的情况%非对角元素HH(ii,jj)=-Uamp(ii)*Uamp(jj)*(G(ii,jj)*sin(theta(ii,jj))-B(ii,jj)*cos(theta(ii,jj)));JJ(ii,jj)=Uamp(ii)*Uamp(jj)*(G(ii,jj)*cos(theta(ii,jj))+B(ii,jj)*sin (theta(ii,jj)));NN(ii,jj)=-Uamp(ii)*(G(ii,jj)*cos(theta(ii,jj))+B(ii,jj)*sin(theta(ii,jj)));LL(ii,jj)=-Uamp(ii)*(G(ii,jj)*sin(theta(ii,jj))-B(ii,jj)*cos(theta(ii,jj)));%对角元素HH(ii,ii)=HH(ii,ii)+Uamp(ii)*Uamp(jj)*(G(ii,jj)*sin(theta(ii,jj))-B(ii,jj)*cos(theta(ii,jj)));JJ(ii,ii)=JJ(ii,ii)-Uamp(ii)*Uamp(jj)*(G(ii,jj)*cos(theta(ii,jj))+B(ii,jj)*sin(theta(ii,jj)) );NN(ii,ii)=NN(ii,ii)-Uamp(jj)*(G(ii,jj)*cos(theta(ii,jj))+B(ii,jj)*sin(theta(ii,jj)));LL(ii,ii)=LL(ii,ii)-Uamp(jj)*(G(ii,jj)*sin(theta(ii,jj))-B(ii,jj)*cos(theta(ii,jj)));endendNN(ii,ii)=NN(ii,ii)-2*Uamp(ii)*G(ii,ii);LL(ii,ii)=LL(ii,ii)+2*Uamp(ii)*B(ii,ii);endpxh(1+2*ng:2:2*(n+ng)-1,1:2:2*n-1)=HH';pxh(1+2*ng:2:2*(n+ng)-1,2:2:2*n)=JJ';pxh(2+2*ng:2:2*(n+ng),1:2:2*n-1)=NN';pxh(2+2*ng:2:2*(n+ng),2:2:2*n)=LL';%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%2、不等式约束的雅克比矩阵%%%%%%%%%%%%%%%%%%%% %g1:电源有功出力上下限约束ag1aP=eye(ng,ng);ag1aQ=zeros(ng,ng);ag1ax=zeros(2*n,ng);%g2:电源无功出力上下限约束ag2aP=zeros(ng,ng);ag2aQ=eye(ng,ng);ag2ax=zeros(2*n,ng);%g3:节点电压幅值上下限约束ag3aP=zeros(ng,n);ag3aQ=zeros(ng,n);ag3ax=zeros(2*n,n);for ii=1:nag3ax(2*ii,ii)=1;end%g4:线路潮流上下限约束ag4aP=zeros(ng,nbr);ag4aQ=zeros(ng,nbr);ag4ax=zeros(2*n,nbr);for ii=1:nfor jj=1:nbrif Nl(jj)==iiag4ax(2*ii-1,jj)=-Uamp(Nl(jj))*Uamp(Nr(jj))*(G(Nl(jj),Nr(jj))*sin(theta(Nl(jj),N r(jj)))-B(Nl(jj),Nr(jj))*cos(theta(Nl(jj),Nr(jj))));ag4ax(2*ii,jj)=Uamp(Nr(jj))*(G(Nl(jj),Nr(jj))*cos(theta(Nl(jj),Nr(jj )))+B(Nl(jj),Nr(jj))*sin(theta(Nl(jj),Nr(jj))))-2*Uamp(Nl(jj))*G(Nl(jj),Nr(jj));endif Nr(jj)==iiag4ax(2*ii-1,jj)=Uamp(Nl(jj))*Uamp(Nr(jj))*(G(Nl(jj),Nr(jj))*sin(theta(Nl(jj),Nr (jj)))-B(Nl(jj),Nr(jj))*cos(theta(Nl(jj),Nr(jj))));ag4ax(2*ii,jj)=Uamp(Nl(jj))*(G(Nl(jj),Nr(jj))*cos(theta(Nl(jj),Nr(jj )))+B(Nl(jj),Nr(jj))*sin(theta(Nl(jj),Nr(jj))));endendendpxg=[ag1aP ag2aP ag3aP ag4aP;ag1aQ ag2aQ ag3aQ ag4aQ;ag1ax ag2ax ag3ax ag4ax];%此即为不等式约束的雅克比矩阵%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%3、对角矩阵%%%%%%%%%%%%%%%%%%%%%%%% L_1Z=zeros(2*ng+n+nbr,2*ng+n+nbr);U_1W=zeros(2*ng+n+nbr,2*ng+n+nbr);for ii=1:2*ng+n+nbrL_1Z(ii,ii)=z(ii)/l(ii);U_1W(ii,ii)=w(ii)/u(ii);end%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%海森伯矩阵%%%%%%%%%%%%%%%%%%%%%%%%%% %将海森伯矩阵分为4块:H1,H2,H3,H4%%%%%%%%%%%%%%%%%%%%%H1%%%%%%%%%%%%%%%%%%%%%%A2=diag(a2);H1=zeros(2*(ng+n),2*(ng+n));H1(1:ng,1:ng)=2*A2;%%%%%%%%%%%%%%%%%%%%H2%%%%%%%%%%%%%%%%%%%%%%H2=zeros(2*(ng+n),2*(ng+n));A=zeros(2*n,2*n);Apb=zeros(2*n,2*n,n);Aqb=zeros(2*n,2*n,n);for ii=1:nfor jj=1:n %元素位置为:1 2if ii~=jj % 3 4%对角线上与ii对应的元素%ApApb(2*ii-1,2*ii-1,ii)=Apb(2*ii-1,2*ii-1,ii)+Uamp(ii)*Uamp(jj)*(G(ii,jj)*cos(theta(ii,jj))+B(ii,jj)*sin(theta(i i,jj)));%对角线处1号元素Apb(2*ii-1,2*ii,ii)=Apb(2*ii-1,2*ii,ii)+Uamp(jj)*(G(ii,jj)*sin(theta(ii,jj))-B(ii,jj)*cos(theta(ii,jj)));%对角线处2号元素%%3号元素与之相等%AqAqb(2*ii-1,2*ii-1,ii)=Aqb(2*ii-1,2*ii-1,ii)+Uamp(ii)*Uamp(jj)*(G(ii,jj)*sin(theta(ii,jj))-B(ii,jj)*cos(theta(ii,jj)));%对角线处1号元素Aqb(2*ii-1,2*ii,ii)=Aqb(2*ii-1,2*ii,ii)-Uamp(jj)*(G(ii,jj)*cos(theta(ii,jj))+B(ii,jj)*sin(theta(ii,jj)));%对角线处2号元素%%3号元素与之相等%对角线上与jj对应的元素%ApApb(2*jj-1,2*jj-1,ii)=Uamp(ii)*Uamp(jj)*(G(ii,jj)*cos(theta(ii,jj))+B(ii,jj)*sin(theta(i i,jj)));%对角线处1号元素Apb(2*jj-1,2*jj,ii)=-Uamp(ii)*(G(ii,jj)*sin(theta(ii,jj))-B(ii,jj)*cos(theta(ii,jj))); %对角线处2号元素Apb(2*jj,2*jj-1,ii)=Apb(2*jj-1,2*jj,ii);%3号元素与2号元素相等%AqAqb(2*jj-1,2*jj-1,ii)=Uamp(ii)*Uamp(jj)*(G(ii,jj)*sin(theta(ii,jj))-B(ii,jj)*cos(theta(ii,jj)));%对角线处1号元素Aqb(2*jj-1,2*jj,ii)=Uamp(ii)*(G(ii,jj)*cos(theta(ii,jj))+B(ii,jj)*sin(theta(ii,jj ))); %对角线处2号元素Aqb(2*jj,2*jj-1,ii)=Aqb(2*jj-1,2*jj,ii);%3号元素与2号元素相等%4号元素为0%非对角线行元素%ApApb(2*ii-1,2*jj-1,ii)=-Uamp(ii)*Uamp(jj)*(G(ii,jj)*cos(theta(ii,jj))+B(ii,jj)*sin(theta(ii,jj)) );%非对角线行处1号元素Apb(2*ii-1,2*jj,ii)=Uamp(ii)*(G(ii,jj)*sin(theta(ii,jj))-B(ii,jj)*cos(theta(ii,jj)));%非对角线行处2号元素Apb(2*ii,2*jj-1,ii)=-Uamp(jj)*(G(ii,jj)*sin(theta(ii,jj))-B(ii,jj)*cos(theta(ii,jj)));%非对角线行处3号元素Apb(2*ii,2*jj,ii)=-(G(ii,jj)*cos(theta(ii,jj))+B(ii,jj)*sin(theta(ii,jj)));%非对角线行处4号元素%AqAqb(2*ii-1,2*jj-1,ii)=-Uamp(ii)*Uamp(jj)*(G(ii,jj)*sin(theta(ii,jj))-B(ii,jj)*cos(theta(ii,jj)));%非对角线行处1号元素Aqb(2*ii-1,2*jj,ii)=-Uamp(ii)*(G(ii,jj)*cos(theta(ii,jj))+B(ii,jj)*sin(theta(ii,jj)));%非对角线行处2号元素Aqb(2*ii,2*jj-1,ii)=Uamp(jj)*(G(ii,jj)*cos(theta(ii,jj))+B(ii,jj)*sin(theta(ii,jj)));%非对角线行处3号元素Aqb(2*ii,2*jj,ii)=-(G(ii,jj)*sin(theta(ii,jj))-B(ii,jj)*cos(theta(ii,jj)));%非对角线行处4号元素%非对角线列元素%ApApb(2*jj-1,2*ii-1,ii)=Apb(2*ii-1,2*jj-1,ii);%非对角线列处1号元素Apb(2*jj-1,2*ii,ii)=Apb(2*ii,2*jj-1,ii);%非对角线列处2号元素Apb(2*jj,2*ii-1,ii)=Apb(2*ii-1,2*jj,ii);%非对角线列处3号元素Apb(2*jj,2*ii,ii)=Apb(2*ii,2*jj,ii);%%非对角线列处4号元素%AqAqb(2*jj-1,2*ii-1,ii)=Aqb(2*ii-1,2*jj-1,ii);%非对角线列处1号元素Aqb(2*jj-1,2*ii,ii)=Aqb(2*ii,2*jj-1,ii);%非对角线列处2号元素Aqb(2*jj,2*ii-1,ii)=Aqb(2*ii-1,2*jj,ii);%非对角线列处3号元素Aqb(2*jj,2*ii,ii)=Aqb(2*ii,2*jj,ii);%%非对角线列处4号元素endend%对角线上与ii对应的元素%ApApb(2*ii,2*ii-1,ii)=Apb(2*ii-1,2*ii,ii);%对角线处3号元素与2号元素相等Apb(2*ii,2*ii,ii)=-2*G(ii,ii);%对角线处4号元素%AqAqb(2*ii,2*ii-1,ii)=Aqb(2*ii-1,2*ii,ii);%对角线处3号元素与2号元素相等Aqb(2*ii,2*ii,ii)=2*B(ii,ii);%对角线处4号元素endfor ii=1:nA=A+Apb(:,:,ii)*y(2*ii-1)+Aqb(:,:,ii)*y(2*ii);endH2(2*ng+1:2*(ng+n),2*ng+1:2*(ng+n))=A;%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%H3%%%%%%%%%%%%%%%%%%%%%%H3=zeros(2*(ng+n),2*(ng+n));A3=zeros(2*n,2*n);Apc=zeros(2*n,2*n,nbr);for ii=1:nbr%对角线上iiApc(2*Nl(ii)-1,2*Nl(ii)-1,ii)=-Uamp(Nl(ii))*Uamp(Nr(ii))*(G(Nl(ii),Nr(ii))*cos(theta(Nl(ii),Nr(ii)))+B( Nl(ii),Nr(ii))*sin(theta(Nl(ii),Nr(ii))));Apc(2*Nl(ii)-1,2*Nl(ii),ii)=-Uamp(Nr(ii))*(G(Nl(ii),Nr(ii))*sin(theta(Nl(ii),Nr(ii)))-B(Nl(ii),Nr(ii))*cos(theta(Nl(ii),Nr(ii))));Apc(2*Nl(ii),2*Nl(ii)-1,ii)=Apc(2*Nl(ii)-1,2*Nl(ii),ii);Apc(2*Nl(ii),2*Nl(ii),ii)=-2*G(Nl(ii),Nr(ii));%对角线上jjApc(2*Nr(ii)-1,2*Nr(ii)-1,ii)=-Uamp(Nl(ii))*Uamp(Nr(ii))*(G(Nl(ii),Nr(ii))*cos(theta(Nl(ii),Nr(ii)))+B( Nl(ii),Nr(ii))*sin(theta(Nl(ii),Nr(ii))));Apc(2*Nr(ii)-1,2*Nr(ii),ii)=Uamp(Nl(ii))*(G(Nl(ii),Nr(ii))*sin(theta(Nl(ii),Nr(ii)))-B(Nl(ii),Nr(ii))*cos(theta(Nl(ii),Nr(ii))));Apc(2*Nr(ii),2*Nr(ii)-1,ii)=Apc(2*Nr(ii)-1,2*Nr(ii),ii);Apc(2*Nr(ii),2*Nr(ii),ii)=0;%非对角线ijApc(2*Nl(ii)-1,2*Nr(ii)-1,ii)=Uamp(Nl(ii))*Uamp(Nr(ii))*(G(Nl(ii),Nr(ii))*cos(theta(Nl(ii),Nr(ii )))+B(Nl(ii),Nr(ii))*sin(theta(Nl(ii),Nr(ii))));Apc(2*Nl(ii)-1,2*Nr(ii),ii)=-Uamp(Nl(ii))*(G(Nl(ii),Nr(ii))*sin(theta(Nl(ii),Nr(ii)))-B(Nl(ii),Nr(ii))*cos(theta(Nl(ii),Nr(ii))));Apc(2*Nl(ii),2*Nr(ii)-1,ii)=Uamp(Nr(ii))*(G(Nl(ii),Nr(ii))*sin(theta(Nl(ii),Nr(ii)))-B(Nl(ii),Nr(ii))*cos(theta(Nl(ii),Nr(ii))));Apc(2*Nl(ii),2*Nr(ii),ii)=G(Nl(ii),Nr(ii))*cos(theta(Nl(ii),Nr(ii))) +B(Nl(ii),Nr(ii))*sin(theta(Nl(ii),Nr(ii)));%非对角线jiApc(2*Nr(ii)-1,2*Nl(ii)-1,ii)=Apc(2*Nl(ii)-1,2*Nr(ii)-1,ii);Apc(2*Nr(ii)-1,2*Nl(ii),ii)=Apc(2*Nl(ii),2*Nr(ii)-1,ii);Apc(2*Nr(ii),2*Nl(ii)-1,ii)=Apc(2*Nl(ii)-1,2*Nr(ii),ii);Apc(2*Nr(ii),2*Nl(ii),ii)=Apc(2*Nl(ii),2*Nr(ii),ii);%求和c=z+w;A3=A3+Apc(:,:,ii)*c(2*ng+n+ii);endH3(2*ng+1:2*(ng+n),2*ng+1:2*(ng+n))=A3;%%%%%%%%%%%%%%%%%%%%%%%H4%%%%%%%%%%%%%%%%%%%%%%%%%H4=pxg*(L_1Z-U_1W)*pxg';%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%H=-H1+H2+H3-H4;%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%形成常数项%%%%%%%%%%%%%%%%%%%%%%%%% %Lyh=zeros(2*n,1);for ii=1:nh(2*ii-1)=Pg(ii)-Pd(ii);h(2*ii)=Qg(ii)-Qd(ii);for jj=1:nh(2*ii-1)=h(2*ii-1)-Uamp(ii)*Uamp(jj)*(G(ii,jj)*cos(theta(ii,jj))+B(ii,jj)*sin(t heta(ii,jj)));h(2*ii)=h(2*ii)-Uamp(ii)*Uamp(jj)*(G(ii,jj)*sin(theta(ii,jj))-B(ii,jj)*cos(theta(ii,jj)));endendLy=h;%Lz%g(x)gx=zeros(2*ng+n+nbr,1);gx(1:ng)=x(1:ng);gx(ng+1:2*ng)=x(ng+1:2*ng);gx(2*ng+1:2*ng+n)=x(2*ng+2:2:2*(ng+n));for ii=1:nbrgx(2*ng+n+ii)=Uamp(Nl(ii))*Uamp(Nr(ii))*(G(Nl(ii),Nr(ii))*cos(theta( Nl(ii),Nr(ii)))+B(Nl(ii),Nr(ii))*sin(theta(Nl(ii),Nr(ii))))-Uamp(Nl(ii))*Uamp(Nl(ii))*G(Nl(ii),Nr(ii));endLz=gx-l-gmin;%LwLw=gx+u-gmax;%Lle=ones(2*ng+n+nbr,1);LZ=zeros(2*ng+n+nbr,2*ng+n+nbr);for ii=1:2*ng+n+nbr;LZ(ii,ii)=l(ii)*z(ii);endLl=LZ*e-miu*e;%LuUW=zeros(2*ng+n+nbr,2*ng+n+nbr);for ii=1:2*ng+n+nbrUW(ii,ii)=u(ii)*w(ii);endLu=UW*e+miu*e;%Lx'Lx1=zeros(2*(ng+n),1);Lx1(1:ng)=2*a2.*x(1:ng)+a1;Lx2=pxh*y;Lx3=pxg*c;Lx41=zeros(2*(ng+n),1);Lx42=zeros(2*(ng+n),1);for ii=1:2*ng+n+nbrLx41(ii)=(Ll(ii)+z(ii)*Lz(ii))/l(ii);Lx42(ii)=(Lu(ii)-w(ii)*Lw(ii))/u(ii);endLx4=pxg*(Lx41+Lx42);Lx=Lx1-Lx2-Lx3;Lxx=Lx+Lx4; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%求出修正量%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %dx,dyHxy=[H pxh;pxh' zeros(2*n,2*n)];LxLy=[Lxx;-Ly];I=eye(2*(ng+n)+2*n);dxdy=I/Hxy*LxLy;dx=dxdy(1:2*(ng+n));dy=dxdy(2*(ng+n)+1:2*(ng+n)+2*n);%dldl=pxg'*dx+Lz;%dudu=-pxg'*dx-Lw;%dzdz=zeros(2*ng+n+nbr,1);for ii=1:2*ng+n+nbrdz(ii)=(-Ll(ii)-z(ii)*dl(ii))/l(ii);end%dwdw=zeros(2*ng+n+nbr,1);for ii=1:2*ng+n+nbrdw(ii)=(-Lu(ii)-w(ii)*du(ii))/u(ii);end %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%计算alfap和alfad%%%%%%%%%%%%%%%%%%%%%%%% alfap=1;alfad=1;for ii=1:2*ng+n+nbrif dl(ii)<0&&-l(ii)/dl(ii)<alfapalfap=-l(ii)/dl(ii);endif du(ii)<0&&-u(ii)/du(ii)<alfapalfap=-u(ii)/du(ii);endif dz(ii)<0&&-z(ii)/dz(ii)<alfadalfad=-z(ii)/dz(ii);endif dw(ii)>0&&-w(ii)/dw(ii)<alfadalfad=-w(ii)/dw(ii);endendalfap=0.9995*alfap;alfad=0.9995*alfad;x=x+alfap*dx;l=l+alfap*dl;u=u+alfap*du;y=y+alfad*dy;z=z+alfad*dz;w=w+alfad*dw;%迭代功率、电压幅值和相角for ii=1:ngPg(Nbus(ii))=x(ii);Qg(Nbus(ii))=x(ng+ii);endfor ii=1:nUamp(ii)=x(2*(ng+ii));Dlta(ii)=x(2*(ng+ii)-1);end %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% k=k+1;elsebreak;endendfcost=0;for ii=1:ngfcost=fcost+a2(ii)*Pg(Nbus(ii))*Pg(Nbus(ii))+a1(ii)*Pg(Nbus(ii))+a0( ii);endfcostkplot(0:k,Gap(1:k+1),':*');PgQgUampfor ii=1:nif Type(ii)==3Dlta=Dlta-Dlta(ii)*ones(n,1);endendDlta。

基于ADMAT自动微分工具箱的内点法最优潮流计算鲍海波1，韦化1,2，莫东1（1．广西大学电力系统最优化研究所，2．广西电力系统最优化与节能技术重点实验室，广西壮族自治区南宁市530004）摘要：为了提高最优潮流程序的通用性和灵活性，避免计算内点法最优潮流时繁琐的导数公式推导，利用ADMAT自动微分工具箱完成了一种更为灵活的内点法最优潮流计算。

无需手动编程构造雅克比矩阵和海森矩阵，减少了编程出错可能性，且便于修改、增减最优潮流模型的目标函数和约束条件。

程序实现过程中，结合ADMAT工具箱的求解特点，研究了约束条件和目标函数的雅克比矩阵和海森矩阵的稀疏模式，对程序进一步优化。

多个系统的测试结果表明，自动微分技术应用于电力系统最优潮流计算的可行性和优越性，所设计的程序具有较高的计算效率。

关键词：自动微分；最优潮流；内点法；ADMAT0 引言经过四十多年的发展，最优潮流（optimal power flow ，OPF）技术在电力系统运行、规划、调度等领域得到了广泛的应用[1-2]。

现今，先后用于最优潮流计算方法有：牛顿法、内点法、遗传算法、差分进化算法等一系列方法[3-9]，其中内点法具有多项式时间复杂性、计算效率高、鲁棒性好等优点，已经成功应用于最优潮流计算的工程实际。

内点法最优潮流计算过程中，需要手动推导目标函数和约束条件的1阶和2阶导数，以获得它们的雅克比（Jacobian）矩阵和海森（Hessian）矩阵。

这种手动编程方式具有很多弊端。

一方面，推导导数计算公式繁琐、容易出错，且代码不易于调试；另一方面，OPF在不同领域应用需要变化模型，增减或者修改部分约束、改变目标函数时，代码改动很大，限制了OPF程序代码的灵活性和可扩展性。

自动微分（automation differentiation，AD）技术的出现，成功解决了这个问题。

AD技术是机械的运用链式求导法则对计算机程序形式的函数求导的一种技术[10]。

基于内点法的最优潮流计算摘要内点法是一种能在可行域内部寻优的方法，即从初始内点出发，沿着中心路径方向在可行域内部直接走向最优解的方法。

其中路径跟踪法是目前最具有发展潜力的一类内点算法，该方法鲁棒性强，对初值的选择不敏感，在目前电力系统优化问题中得到了广泛的应用。

本文采用路径跟踪法进行最优求解，首先介绍了路径跟踪法的基本模型，并且结合具体算例，用编写的Matlab程序进行仿真分析，验证了该方法在最优潮流计算中的优越性能。

关键词：最优潮流、内点法、路径跟踪法、仿真目次0、引言 (1)1、路径跟踪法的基本数学模型 (2)2、路径跟踪法的最优潮流求解思路 (3)3、具体算例及程序实现流程 (7)3.1、算例描述 (7)3.2、程序具体实现流程 (8)4、运行结果及分析 (12)4．1 运行结果 (12)4．2结果分析 (18)5、结论 (19)6、编程中遇到的问题 (20)参考文献 (22)附录 (23)0、引言电力系统最优潮流，简称OPF（Optimal Power Flow）。

OPF问题是一个复杂的非线性规划问题，要求满足待定的电力系统运行和安全约束条件下，通过调整系统中可利用控制手段实现预定目标最优的系统稳定运行状态。

针对不同的应用，OPF模型课以选择不同的控制变量、状态变量集合，不同的目标函数，以及不同的约束条件，其数学模型可描述为确定一组最优控制变量u，以使目标函数取极小值，并且满足如下等式和不等式。

（0-1）其中为优化的目标函数，可以表示系统运行成本最小、或者系统运行网损最小。

为等式约束，表示满足系统稳定运行的功率平衡。

为不等式约束，表示电源有功出力的上下界约束、节点电压上下线约束、线路传输功率上下线约束等等。

电力系统最优潮流算法大致可以分为两类：经典算法和智能算法。

其中经典算法主要是指以简化梯度法、牛顿法、内点法和解耦法为代表的基于线性规划和非线性规划以及解耦原则的算法，是研究最多的最优潮流算法, 这类算法的特点是以一阶或二阶梯度作为寻找最优解的主要信息。

基于内点法的电力系统最优潮流算法研究毕业论文目录摘要....................................... 错误！未定义书签。

ABSTRACT ....................................... 错误！未定义书签。

目录.. (1)1绪论 (3)1.1引言 (3)1.2电力系统最优潮流计算的发展历史及现状 (3)1.3本文所做工作 (4)2、电力系统最优潮流算法介绍 (6)2.1最优潮流计算的基本数学模型 (6)2.1.1目标函数 (6)2.1.2等式约束条件 (7)2.1.3不等式约束条件 (7)2.2电力系统最优潮流的算法简介 (8)2.2.1 线性规划法 (8)2.2.2 二次规划法 (8)2.2.3 牛顿法 (9)2.2.4内点法 (9)2.2.5电力系统最优潮流计算的新兴算法 (10)3、原对偶内点法 (11)3.1原对偶内点法的数学原理 (11)3.2目标函数的收敛条件 (14)3.3 初值的选取 (15)3.4 利用原对偶内点法进行潮流计算的方法 (15)4.基于原对偶内点法的电力系统最优潮流计算 (17)4.1电力系统最优潮流计算中的各项数学模型 (17)4.1.1最优潮流计算的目标函数 (17)4.1.2最优潮流计算的等式约束条件 (17)4.1.3最优潮流计算的不等式约束条件 (18)4.2各项数学模型的具体表达 (18)4.2.1目标函数的各偏导数及相应矩阵 (18)4.2.2等式约束的各偏导数及相应矩阵 (19)4.2.3不等式约束的各偏导数及相应矩阵 (23)4.2.4对模型中各节点的不等式约束条件的处理 (26)4.3算例分析 (27)4.3.1 MATLAB简介 (27)4.3.2具体的计算流程 (27)4.3.3 IEEE-14标准测试系统运算结果 (28)5总结与展望 (31)5.1本文总结 (31)5.2今后展望 (31)参考文献 (32)附录IEEE-14标准测试系统数据 (34)致谢....................................... 错误！未定义书签。

内点法在具有新颖数据结构的最优潮流问题中的应用人名？摘要：基于原问题的扰动KKT条件，本文提出了一种新的内点非线性规划算法（OPF），以求解最优潮流问题。

通过中心路径方向概念，我们将该算法扩展用于普通潮流（PF）和近似最优潮流问题。

对于后者，CPU时间，可大幅度减少。

为了有效处理函数不等式约束，推导了简约修正方程，该方程的结构取决于等式约束的结构。

通过重新排列修正方程，提出了一种新颖的数据结构。

与牛顿法OPF的传统数据结构相比，提出方案的填充元大约减少了一半，CPU时间减少约15％。

该算法包括四个目标函数和两种不同的数据结构。

在14到1047节点系统上，广泛的数字仿真表明：由于提出方法的鲁棒性和快速的执行时间，非常有前景应用于大规模系统。

关键词：最优潮流，内点非线性规划，扰动KKT条件，中心方向，近似最优潮流一、引言以下自己对着翻译最优潮流(潮流) 由卡尔庞捷[1] 在早期的20 世纪60 年代中定义的因为它一直吸引着许多研究者作为一种潜在的强大的工具为电力系统运营商和规划。

过去三年，各种的优化方法，如线性和非线性的二次编程（LP、QP 和NP），分解和牛顿法[2-6]，已应用于此吸引人的话题。

但是，虽然电力系统正变得越来越大、越来越复杂，潮流问题不过正变得越来越困难。

到目前为止，许多开放问题尚待解决。

求解LP 问题的多项式时间内点法(IPM) 引起了研究人员的关注。

此外解决了在电力系统中的各种优化问题的IPM [7-10]。

它是向IPM QP 和NP问题的必然趋势。

非线性申请[11-16] 报告了有价值的结果。

虽然IPM 非凸NP问题的理论研究已被发展，与大型应用程序关联的许多问题未被认识清楚不过。

从理论上讲，IPM 算法可以快速聚合到一个可行，但不是最佳的解决方案，沿中心方向[18] （请参见11 B 部分）。

虽然这个方向不能确定最佳的解决方案，其快速收敛是很有吸引力。

无法识别的问题是居中的方向可以提供任何有价值的应用程序，PF 和0PF-此外，优秀的IPM 算法必须具备适合潮流的数据结构，因为它强烈影响算法的性能。

内点最优潮流算法自动微分的有效执行技术摘要：本文提出了一种改进的内点矩形最优潮流（OPF)算法的自动微分（AD）技术执行过程。

有别于现有的AD技术执行过程，该执行过程增加了一个识别由AD技术生成的所有定常一阶和二阶导数的子程序，并在迭代前生成一个定常导数列表。

在内点OPF算法的每次迭代工程中，只通过AD工具更新变化的导数。

ADC这一优秀的软件作为AD的一个基本工具，完成上述执行工程。

AD技术结合用户自定义模型界面，增强了计算性能和灵活性。

大规模的电力系统算例研究表明，该算法在保持代码可维护性、灵活性的同时,计算速度接近手动编程。

这篇文章证明，AD技术具有应用于电力系统在线运行环境的潜力，可取代传统手动编程求导，大大减轻软件开发人员的负担。

关键字：ADC，自动微分，内点法，操作符重载，最优潮流1.引言近年来，内点法（IPM)凭借其出色的计算性能和鲁棒性，已经被广泛应用于求解大规模电力系统最优潮流（OPF）问题。

在内点法OPF 中,计算目标函数与约束条件的梯度，雅可比(J acobian) 矩阵和海森( Hessian) 矩阵是很重要的部分。

为了获得上述矩阵,开发者不得不手动推导一阶和二阶导数计算公式并手动编程。

这种手动编程方式具有以下缺点: ①推导导数计算公式过于繁琐且易于出错; ②将上述公式手动编程并调试工作量大且容易出错; ③当加入新设备或复杂装置(如柔性交流输电系统( FACTS) 和高压直流( HVDC) 装置)时，增减或修改约束条件、改变目标函数时会很繁琐。

自动微分(AD) 技术的使用克服了手动编程的缺点,与其他微分方法(如数值差分、符号微分) 相比,AD 避免了截断误差,对中央处理器(CPU) 时间和内存空间的占用都远小于上述方法。

文献5中在电力系统动态仿真中采用AD计算jacobian矩阵。

文献6-8采用AD算法计算电力系统潮流。

文献9在计算电力系连续潮流jacobian矩阵和灵敏度时采用AD。