第六章 信頼域方法

Posted 2021-06-29 是璇子鸭

内容来自马昌凤编著的《最优化方法及其Matlab程序设计》，文章仅为个人的学习笔记，感兴趣的朋友详见原书。

信頼域方法和线搜索技术都是在优化算法中求出每次迭代的位移，从而确定新的迭代点。不同的是：线搜索技术是先产生搜索方向，然后确定搜索步长；而信頼域技术则是直接确定位移，产生新的迭代点。

基本思想：首先给定一个“信頼域半径”作为位移长度的上界，并以当前迭代点为中心，此“上界”为半径“画地为牢”，确定一个被称为“信頼域”的闭球区域。再以此为基础，一步步求解、缩放半径，直至满足迭代终止条件。

信頼域方法的基本结构

信頼域子问题的求解

程序：利用光滑牛顿法求解信頼域子问题，一般适用于（近似）Hess矩阵正定的情况

function  [d,val,lam,k]=trustq(gk,Bk,dta)
% 功能: 求解信赖域子问题:  min qk(d)=gk'*d+0.5*d'*Bk*d, s.t. ||d||<=delta
%输入:  gk是xk处的梯度, Bk是第k次近似Hesse阵, dta是当前信赖域半径
%输出:  d, val分别是子问题的最优点和最优值, lam是乘子值, k是迭代次数.
n=length(gk);   gamma=0.05; 
epsilon=1.0e-6;  rho=0.6;  sigma=0.2;   
mu0=0.05;  lam0=0.05;
d0=ones(n,1);   z0=[mu0,lam0,d0']'; 
u0=[mu0,zeros(1,n+1)]';
k=0;   %k为迭代次数
z=z0;
mu=mu0; lam=lam0; d=d0;
while ( k<=150)     %Step1 of the algorithm
    dh=dah(mu,lam,d,gk,Bk,dta);
     if(norm(dh)<epsilon)
         break;
     end
     A=JacobiH(mu,lam,d,Bk,dta); b=beta(mu,lam,d,gk,Bk,dta,gamma)*u0-dh;
     B=inv(A);      dz=B*b;
     dmu=dz(1); dlam=dz(2); dd=dz(3:n+2);
     m=0;  mk=0;
     while (m<20)
         dhnew=dah(mu+rho^m*dmu,lam+rho^m*dlam,d+rho^m*dd,gk,Bk,dta);
         if(norm(dhnew)<=(1-sigma*(1-gamma*mu0)*rho^m)*dh)
             mk=m;           
             break;
         end
         m=m+1;
     end
     alpha=rho^mk;
     mu=mu+alpha*dmu;
     lam=lam+alpha*dlam;
     d=d+alpha*dd;
     k=k+1;
end
 val=gk'*d+0.5*d'*Bk*d;
%%%%%%%%%%%%%%%%%%%%%%%%%%
function p=phi(mu,a,b)
p=a+b-sqrt((a-b)^2+4*mu);
%%%%%%%%%%%%%%%%%%%%%%%%%%
function dh=dah(mu,lam,d,gk,Bk,dta)
n=length(d);
dh(1)=mu;  dh(2)=phi(mu,lam, dta^2-norm(d)^2);
mh=(Bk+lam*eye(n))*d+gk;
for(i=1:n)
    dh(2+i)=mh(i);
end
dh=dh(:);
%%%%%%%%%%%%%%%%%%%%%%%%%%
function bet=beta(mu,lam,d,gk,Bk,dta,gamma)
dh=dah(mu,lam,d,gk,Bk,dta);
bet=gamma*norm(dh)*min(1,norm(dh));
%%%%%%%%%%%%%%%%%%%%%%%%%%
function A=JacobiH(mu,lam,d,Bk,dta)
n=length(d);
A=zeros(n+2,n+2);
pmu=-4*mu/sqrt((lam+norm(d)^2-dta^2)^2+4*mu^2);
thetak=(lam+norm(d)^2-dta^2)/sqrt((lam+norm(d)^2-dta^2)^2+4*mu^2);
A=[1,0,zeros(1,n);pmu,1-thetak,-2*(1+thetak)*d';zeros(n,1),  d,Bk+lam*eye(n)];

信頼域方法的Matlab程序

程序：牛顿型信頼域方法，其中涉及到对子程序trustq的调用（见上方）

function [xk,val,k]=trustm(x0)
%功能: 牛顿型信赖域方法求解无约束优化问题 min f(x)
%输入: x0是初始迭代点
%输出: xk是近似极小点, val是近似极小值, k是迭代次数
n=length(x0);  x=x0; dta=1;
eta1=0.15; eta2=0.75;  dtabar=2.0;
tau1=0.5; tau2=2.0; epsilon=1e-6;
k=0;  Bk=Hess(x);  %Bk=eye(n);  
while(k<150)
    gk=gfun(x);   
    if(norm(gk)<epsilon)
        break;
    end
    [d,val,lam,ik]=trustq(gk,Bk,dta);
    deltaq=-qk(x,d);
    deltaf=fun(x)-fun(x+d);
    rk=deltaf/deltaq;
    if(rk<=eta1)
        dta=tau1*dta;
    else if (rk>=eta2&norm(d)==dta)
            dta=min(tau2*dta,dtabar);
        else
            dta=dta;
        end
    end
    if(rk>eta1)
        x0=x;     x=x+d;    
       % sk=x-x0;  yk=gfun(x)-gfun(x0);  
        %vk=sqrt(yk'*Bk*yk)*(sk/(sk'*yk)-Bk*yk/(yk'*Bk*yk));
        %Bk=Bk-Bk*yk*yk'*Bk/(yk'*Bk*yk)+sk*sk'/(sk'*yk)+vk*vk'
        %pause
        Bk=Hess(x);
    end
    k=k+1;
end
xk=x;
val=fun(xk);