FR共轭梯度法实验及其程序代码

实验报告2FR共轭梯度法1实验目的掌握外法函数法2实验内容外罚函数法3算法设计①编写主函数：function[min,answ]=OuterPenaltyFunction(x0,eps)实现外罚函数法②定义functiony=PFun(x)用以实现求目标函数的一阶导数；③定义测试函数:functiony=Fun(x)。4程序代码（一）外罚函数法程序代码function[min,answ]=OuterPenaltyFunction(x0,eps)%%----该程序用以实现外罚函数法%%Input：%x0----初始点%eps----搜索精度%%output：%min----目标函数最优值%answ----最优搜索点symsx1x2xk=x0;c=10;%罚因子放大系数33deta=0.1;%罚因子while(1)fun=Fun([x1x2])+deta*PFun([x1,x2]);%惩罚函数f=inline(fun);[min,answ]=FR_Conjugate_gradient(f,xk,eps);%调用FR迭代法函数求确定罚函数最优解ifdouble(deta*PFun(answ))epsbreakenddeta=c*deta;xk=answ;endansw=xk;min=double(Fun(xk));functiony=PFun(x)%%----该函数用以定义惩罚项函数%%Input：%x----自变量%%output：%y----因变量symsx1x2symsx1x2x1=x(1);%x2=x(2);y=(x1+1)^2;functiony=Fun(x)%%----该函数用以定义测试函数%%Input：%x----自变量%%output：%y----因变量symsx1x2x1=x(1);x2=x(2);y=x1^2+x2^2;（二）FR共轭梯度法程序代码：function[min,answ]=FR_Conjugate_gradient(Fun,x0,eps)%%----改程序用以实现FR共轭梯度法%%Input：%x0----初始点%eps----搜索精度%%output：%min----目标函数最优值%answ----最优搜索点symsbetaalphk=1;%迭代次数xk=x0;temp=1;while(norm(DFun(Fun,xk))eps)if(k==1)pk=-DFun(Fun,xk);elsebeta=DFun(Fun,xk)'*DFun(Fun,xk)/temp;pk=-DFun(Fun,xk)+beta*pk;endtemp=DFun(Fun,xk)'*DFun(Fun,xk);temp2=xk+alph*pk';fun=Fun(temp2(1),temp2(2));fun=inline(fun);[min,answ]=OneDimentionSearch(fun,0,0.08,eps);%调用精确一维搜索求步长xk=xk+answ*pk';%answ即为所求的步长k=k+1;endansw=xk;min=double(Fun(xk(1),xk(2)));functiony=DFun(Fun,x)%%----改函数用以实现求目标函数的一阶导数%%Input：%x----求导点%%output：%y----目标函数的一阶导数在x点的值symsx1x2f=Fun(x1,x2);f1=diff(f,x1);f2=diff(f,x2);f1=subs(f1,{x1,x2},x);f2=subs(f2,{x1,x2},x);y=[double(f1),double(f2)]';（三）精确一维搜索程序代码：function[min,answ]=OneDimentionSearch(Fun,x0,step,eps)%%----改程序用以实现一维搜索算法%%Input：%Fun----一维精确搜索函数%x0----初始点%step----搜索步长%eps----搜索精度%%output：%min----目标函数最优值%answ----最优搜索点%%%求初始区间[a,b]=AdvanceAndRetreat(Fun,x0,step);%黄金分割法实现一维搜索[min,answ]=GoldSection(Fun,a,b,eps);function[min,answ]=GoldSection(Fun,a,b,eps)%%----改函数用以实现黄金搜索法%%Input：%Fun----一维精确搜索函数%a----搜索区间的左端点%b----搜索间隔的右端点%eps----搜索精度%%output：%min----目标函数最优值%wnsw----最优搜索点%%formatlongx1=a+.382*(b-a);x2=a+.618*(b-a);while(abs(b-a)eps)f1=Fun(x1);f2=Fun(x2);%比较判断两个分割点处的函数值，进而缩短区间长度if(f1f2)a=x1;x1=x2;x2=a+.618*(b-a);elseif(f1==f2)a=x1;b=x2;x1=a+.382*(b-a);x2=a+.618*(b-a);elseb=x2;x2=x1;x1=a+.382*(b-a);endend%返回搜索点和搜索值answ=(a+b)/2;min=Fun(answ);function[a,b]=AdvanceAndRetreat(Fun,x0,step)%%----改函数用以实现进退法求初始区间%%Input：%Fun----一维精确搜索函数%x0----初始点%step----搜索间隔%%output：%a----搜索区间的左端点%b----搜索间隔的右端点x1=x0+step;ifFun(x1)=Fun(x0)while(1)step=2*step;x2=x1+step;ifFun(x1)=Fun(x2)break;elsex0=x1;x1=x2;endenda=x0;b=x2;elsewhile(1)step=2*step;x2=x0-step;ifFun(x0)=Fun(x2)break;elsex1=x0;x0=x2;endenda=x2;b=x1;end5运行结果测试函数为：01..)(min12221xtsxxxf该函数是教材第147页例4.2.2，其结果与书上的结果近似相等（因为存在微小迭代误差）。理论最优解为：Tx)0,1(*，最优值为:1*f。运行结果如下图：6结果分析该程序实现了外罚函数法解决最优化问题，但是具有局限性，只能求解含有两个自变量的优化问题，原因是在编写FR共轭梯度法程序时的求导函数时人为地将在变量限制在两个。然而，程序代码编写过程思路清晰，可读性强，只要稍微将FR共轭梯度法中的求导函数DFun（x）稍稍修改便可以实现所有功能，即程序具有很强的扩展性。以上运行结果除微小误差外完全与理论解吻合，由此可见程序的正确性。