优化求解基于matlab内点法求解实时电价最优问题含Matlab源码 1161期

Posted 紫极神光

tags:

篇首语:本文由小常识网(cha138.com)小编为大家整理,主要介绍了优化求解基于matlab内点法求解实时电价最优问题含Matlab源码 1161期相关的知识,希望对你有一定的参考价值。

一、简介


法。而且无论是面对LP还是QP,内点法都显示出了相当的极好的性能,例如多项式的算法复杂度。

二、源代码

clear;
clc;
errArr=[];
%%
%初始化!!!
initial;
% Start clock
t1 = clock;
%%
ROU=sl'*MU_MIN+su'*MU_MAX;
MUt=SIGMA*ROU/(2*length(sl));%初始对偶因子与惩罚因子计算%
ik=0;%计迭代次数!!!
%迭代循环过程!!
while(abs(ROU)>=err)
    %%
    %Calcute h,g matrix
    ROU=sl'*MU_MIN+su'*MU_MAX;
        errArr=[errArr;ROU;];
    SIGMA=0;
    MU=SIGMA*ROU/(2*length(sl));         %中心参数置零%
    for i=1:30
        temp=0;
        for j=1:30
            temp=temp-V(j)*aY(i,j)*cos(Vth(i)-Vth(j)-Yth(i,j));
        end
        if (i>6) 
            tPg=0;
        else
            tPg=Pg(i);
        end
        h(i)=tPg-Pd(i)+V(i)*temp;
    end
    for i=1:30
        temp=0;
        for j=1:30
            temp=temp-V(j)*aY(i,j)*sin(Vth(i)-Vth(j)-Yth(i,j));
        end
        if (i>6) 
            tQg=0;
        else
            tQg=Qg(i);
        end
        h(i+30)=tQg-Qd(i)+V(i)*temp;
    end        % Cal  h END
    
    for i=1:6
        g(i)=Pg(i);
        g(i+6)=Qg(i);
    end
    for i=1:30
        g(i+12)=V(i);
    end        % Cal  g  END
    %Calcute h,g matrix END
    %%
    %Calculate Jacobian&Hessian matix
    %First Step: Jf,Hf
    for i=1:6
        Jf(i)=2*gencost(i,5)*Pg(i)+gencost(i,6);
        Hf(i,i)=2*gencost(i,6);
    end
    %Second Step: Jh, h为等式约束
    for i=1:6         %6行对Pg求导,由此已求出
        Jh(i,i)=1; 
    end
    for i=7:12        %7-12行对Qg求导,由此已求出
        Jh(i,i+24)=1; 
    end
    for i=1:30       %形成13-42行的1-60for j=1:30
            tempVp=0;
            tempVq=0;
            if (j==i)
                for k=1:30
                    tempVp=tempVp-V(k)*aY(j,k)*cos(Vth(j)-Vth(k)-Yth(j,k));
                    tempVq=tempVq-V(k)*aY(j,k)*sin(Vth(j)-Vth(k)-Yth(j,k));
                end
                Jh(12+j,i)=tempVp-aY(j,j)*V(j)*cos(Yth(j,j));
                Jh(12+j,30+i)=tempVq+aY(j,j)*V(j)*sin(Yth(j,j));
            else
                Jh(12+j,i)=-aY(i,j)*V(i)*cos(Vth(i)-Vth(j)-Yth(i,j));
                Jh(12+j,30+i)=-aY(i,j)*V(i)*sin(Vth(i)-Vth(j)-Yth(i,j));
            end
        end
    end
    for i=1:30       %形成43-72行的1-60for j=1:30
            tempVp=0;
            tempVq=0;
            if (j==i)
                for k=1:30
                    tempVp=tempVp+aY(j,k)*V(k)*sin(Vth(j)-Vth(k)-Yth(j,k));
                    tempVq=tempVq-aY(j,k)*V(k)*cos(Vth(j)-Vth(k)-Yth(j,k));
                end
                tempVp=tempVp-V(j)*aY(j,j)*sin(-Yth(j,j));
                tempVq=tempVq+V(j)*aY(j,j)*cos(-Yth(j,j));
                Jh(42+j,i)=V(i)*tempVp;
                Jh(42+j,30+i)=V(i)*tempVq;
            else
                Jh(42+j,i)=-aY(i,j)*V(i)*V(j)*sin(Vth(i)-Vth(j)-Yth(i,j));
                Jh(42+j,30+i)=aY(i,j)*V(i)*V(j)*cos(Vth(i)-Vth(j)-Yth(i,j));
            end
        end
    end
    %Third Step: Hh
    %有功部分
    for i=1:30
        for j=1:30
            for k=j:30
                if (j==k&&i~=j) 
                    Hh(j+12,k+12,i)=0;          %VV
                    Hh(j+42,k+42,i)=V(i)*aY(i,j)*V(j)*cos(Vth(i)-Vth(j)-Yth(i,j)); %%thth
                elseif (j==k&&i==j)
                    Hh(j+12,k+12,i)=-2*aY(j,j)*cos(Yth(i,i));  %VV
                    temp=0;                   %thth
                    for l=1:30
                        temp=temp+aY(j,l)*V(l)*cos(Vth(j)-Vth(l)-Yth(j,l));
                    end
                    temp=temp-aY(i,i)*V(i)*cos(-Yth(i,i));
                    Hh(j+42,k+42,i)=V(i)*temp;
                elseif (k==i)
                    Hh(j+12,k+12,i)=-aY(i,j)*cos(Vth(i)-Vth(j)-Yth(i,j));    %VV
                    Hh(k+12,j+12,i)=Hh(j+12,k+12,i);
                    Hh(j+42,k+42,i)=-V(i)*aY(i,j)*V(j)*cos(Vth(i)-Vth(j)-Yth(i,j));  %thth
                    Hh(k+42,j+42,i)=Hh(j+42,k+42,i);
                elseif (j==i)
                    Hh(j+12,k+12,i)=-aY(i,k)*cos(Vth(i)-Vth(k)-Yth(i,k));    %VV
                    Hh(k+12,j+12,i)=Hh(j+12,k+12,i);
                    Hh(j+42,k+42,i)=-V(i)*aY(i,k)*V(k)*cos(Vth(i)-Vth(k)-Yth(i,k));  %thth
                    Hh(k+42,j+42,i)=Hh(j+42,k+42,i);
                end
            end
        end
    end                 %至此已形成(13-4213-42)和(42-7243-72for i=1:30
        for j=1:30
            for k=1:30
                if (j==k&&i~=j) 
                    Hh(j+42,k+12,i)=-V(i)*aY(i,j)*sin(Vth(i)-Vth(j)-Yth(i,j));      %thV
                elseif (j==k&&i==j)
                    temp=0;                   %thV
                    for l=1:30
                        temp=temp+aY(j,l)*V(l)*sin(Vth(j)-Vth(l)-Yth(j,l));    
                    end
                    Hh(j+42,k+12,i)=temp-V(i)*aY(i,i)*sin(-Yth(i,i));
                elseif (j==i)
                    Hh(j+42,k+12,i)=V(i)*aY(i,k)*sin(Vth(i)-Vth(k)-Yth(i,k));   %thV   
                elseif (k==i)
                    Hh(j+42,k+12,i)=-V(j)*aY(i,j)*sin(Vth(i)-Vth(j)-Yth(i,j));  %thV
                end
            end
        end
        Hh(13:42,43:72,i)=Hh(43:72,13:42,i)';
    end                           %至此已形成(42-7213-42)和(13-4243-72%无功部分
    for i=1:30
        for j=1:30
            for k=j:30
                if (j==k&&i~=j) 
                    Hh(j+12,k+12,i+30)=0;          %VV
                    Hh(j+42,k+42,i+30)=V(i)*aY(i,j)*V(j)*sin(Vth(i)-Vth(j)-Yth(i,j)); %%thth
                elseif (j==k&&i==j)
                    Hh(j+12,k+12,i+30)=2*aY(j,j)*sin(Yth(i,i));  %VV
                    temp=0;                   %thth
                    for l=1:30
                        temp以上是关于优化求解基于matlab内点法求解实时电价最优问题含Matlab源码 1161期的主要内容,如果未能解决你的问题,请参考以下文章

C++内点法求解大规模线性规划问题——对标MATLAB中linprog函数

C++内点法求解大规模线性规划问题——对标MATLAB中linprog函数

C++内点法求解大规模线性规划问题——对标MATLAB中linprog函数

C++内点法求解大规模线性规划问题——对标MATLAB中linprog函数

优化求解基于matlab遗传算法求解峰谷分时电价电动汽车充电负荷优化问题含Matlab源码 210期

第九章 罚函数法