文档介绍：function[x,mu,lambda,output]=multphr(fun,hf,gf,dfun,dhf,dgf,x0)%功能:用乘子法解一般约束问题:minf(x),(x)=0,g(x).=0%输入:x0是初始点,fun,dfun分别是目标函数及其梯度;%hf,dhf分别是等式约束(向量)函数及其Jacobi矩阵的转置;%gf,dgf分别是不等式约束(向量)函数及其Jacobi矩阵的转置;%输出:x是近似最优点,mu,lambda分别是相应于等式约束和不%等式约束的乘子向量;output是结构变量,输出近似极小值f,迭%代次数,内迭代次数等maxk=500;%最大迭代次数sigma=;%罚因子eta=;theta=;%PHR算法中的实参数k=0;ink=0;%k,ink分别是外迭代和内迭代次数epsilon=1e-5;%终止误差值x=x0;he=feval(hf,x);gi=feval(gf,x);n=length(x);l=length(he);m=length(gi);%选取乘子向量的初始值mu=*ones(l,1);lambda=*ones(m,1);btak=10;btaold=10;%用来检验终止条件的两个值while(btak>epsilon&k<maxk)%调用BFGS算法程序求解无约束子问题[x,ival,ik]=bfgs('mpsi','dmpsi',x0,fun,hf,gf,dfun,dhf,dgf,mu,lambda,sigma);ink=ink+ik;he=feval(hf,x);gi=feval(gf,x);btak=;for(i=1:l),btak=btak+he(y2;endfor(i=1:m)temp=min(gi(i),lambda(i)/sigma);btak=btak+tempA2;endbtak=sqrt(btak);ifbtak>epsilonif(k>=2&btak>theta*btaold)sigma=eta*sigma;end%更新乘子向量for(i=1:l),mu(i)=mu(i)-sigma*he(i);endfor(i=1:m)lambda(i)=max(,lambda(i)-sigma*gi(i));endendk=k+1;btaold=btak;x0=x;endf=feval(fun,x);=f;=k;=ink;=btak;%%%%%%%%%%%%%%%%%%增广拉格朗日函数%%%%%%%%%%%%%%%%%%%%%functionpsi=mpsi(x,fun,hf,gf,dfun,dhf,dgf,mu,lambda,sigma)f=feval(fun,x);he=feval(hf,x);gi=feval(gf,x);l=length(he);m=length(gi);psi=f;s1=;for(i=1:l)psi=psi-he(i)*mu(i);s1=s1+he(i)A2;endpsi=psi+*sigma*s1;s2=;for(i=1:m)s3=max(,lambda(i)-sigma*gi(i));s2=s2