Matlab最優(yōu)化編程例子

1、題目：分別用最速下降法、FR共軛梯度法、DFP法和BFGS法求解問題:minf(x) = X2 2x x + 4x2 + x 一3x1 1 2 2 1 2取初始點(diǎn)X二(1,1，通過Matlab編程實(shí)現(xiàn)求解過程。公用函數(shù)如下：1、function f= fun( X )%所求問題目標(biāo)函數(shù)f=X(1)入2-2*X(1)*X(2)+4*X(2)入2+X(1)-3*X(2);end2、function g= gfun( X )%所求問題目標(biāo)函數(shù)梯度 g=2*X(1)-2*X(2)+1,-2*X(1)+8*X(2)-3;end3、function He = Hess( X )%所求問題目標(biāo)函數(shù)Hesse

2、矩陣 n=length(X);He=zeros(n,n);He=2,-2;-2,4;End解法一：最速下降法function x,val,k = grad( fun,gfun,x0 )%功能：用最速下降法求無約束問題最小值%輸入：xO是初始點(diǎn)，fun和gfun分別是目標(biāo)函數(shù)和梯度%輸出：x、val分別是最優(yōu)點(diǎn)和最優(yōu)值，k是迭代次數(shù)maxk=5000;%最大迭代次數(shù)rho=O.5;sigma=O.4;k=O;eps=1Oe-6;while(kmaxk)g=feval(gfun,xO) %計(jì)算梯度d=-g %計(jì)算搜索方向if(norm(d)eps)break;endm=O;mk=O;while(m

3、2O)if(feval(fun,x0+rhoAm*d)feval(fun,x0)+sigma*rhoAm*g*d) mk=break;endm=m+1;endx0=x0+rhoAmk*d;k=k+1;endx=x0;val=feval(fun,x0);end解法二：FR共軛梯度法function x,val,k = frcg( fun,gfun,x0 )%功能：用FR共軛梯度法求無約束問題最小值%輸入：xO是初始點(diǎn)，fun和gfun分別是目標(biāo)函數(shù)和梯度 %輸出：x、val分別是最優(yōu)點(diǎn)和最優(yōu)值，k是迭代次數(shù) maxk=5O O O ; %最大迭代次數(shù)rho=O.5;sigma=O.4;k=O;e

4、ps=1Oe-6;n=length(xO);while(k=O.O)d=-g;endendif(norm(g)eps)break;endm=O;mk=O;while(m20)if(feval(fun,x0+rhoAm*d)feval(fun,x0)+sigma*rhoAm*g*d) mk=break;endm=m+1;endxO=xO+rhoFk*d;val=feval(fun,x0);g0=g;d0=d;k=k+1;endx=x0;val=feval(fun,x0);end解法三：DFP法function x,val,k = dfp( fun,gfun,x0 )%功能：用DFP法求無約束問題

5、最小值%輸入：xO是初始點(diǎn)，fun和gfun分別是目標(biāo)函數(shù)和梯度 %輸出：x、val分別是最優(yōu)點(diǎn)和最優(yōu)值，k是迭代次數(shù) maxk=5000;%最大迭代次數(shù)rho=O.5;sigma=O.4;k=O;eps=1Oe-6;n=length(xO);Hk=inv(feval(Hess,xO); while(kmaxk)gk=feval(gfun,xO);if(norm(gk)eps)break;enddk=-Hk*gk;dk=dk;m=O;mk=O;while(m2O)if(feval(fun,x0+rhoAm*dk)0) Hk=Hk-(Hk*yk)*yk)*Hk)/(yk*Hk*yk)+(sk*s

6、k)/(sk*yk);endk=k+1;x0=x;endval=feval(fun,x0);end解法四：BFGS法function x,val,k = bfgs( fun,gfun,x0 ) %功能：用BFGS法求無約束問題最小值%輸入：xO是初始點(diǎn)，fun和gfun分別是目標(biāo)函數(shù)和梯度 %輸出：x、val分別是最優(yōu)點(diǎn)和最優(yōu)值，k是迭代次數(shù) maxk=5000;%最大迭代次數(shù)rho=O.5;sigma=O.4;k=O;eps=1Oe-6;n=length(xO);Bk=eye(n);while(kmaxk)gk=feval(gfun,xO);if(norm(gk)eps)break;enddk=-Bk*gk;m=0;mk=0;while(m20)new=sigma*rhoF*gk*dk;old=feval(fun,x0);if (feval(fun,x0+