stepan
/
dynare
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Slightly simplified checks in computing gradient and Hessian for optimizer = 5.

time-shift
Marco Ratto 2011-02-14 11:51:15 +01:00
parent 13ea1c0046
commit 7fb471e9cc
3 changed files with 47 additions and 34 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

14
matlab/mr_gstep.m
View File

@ -27,7 +27,8 @@ n=size(x,1);
if init, if init,
gstep_ = options_.gstep; gstep_ = options_.gstep;
h1=max(abs(x),sqrt(gstep_)*ones(n,1))*eps^(1/4); % h1=max(abs(x),sqrt(gstep_)*ones(n,1))*eps^(1/4);
h1=options_.gradient_epsilon*ones(n,1);
return return
end end
if nargin<4, if nargin<4,
@ -71,7 +72,8 @@ while i<n,
icount = 0; icount = 0;
h0=h1(i); h0=h1(i);
while (abs(dx(it))<0.5*htol | abs(dx(it))>(2*htol)) & icount<10 & ic==0, % while (abs(dx(it))<0.5*htol | abs(dx(it))>(2*htol)) & icount<10 & ic==0,
while (abs(dx(it))<0.5*htol) & icount<10 & ic==0,
%while abs(dx(it))<0.5*htol & icount< 10 & ic==0, %while abs(dx(it))<0.5*htol & icount< 10 & ic==0,
icount=icount+1; icount=icount+1;
if abs(dx(it)) ~= 0, if abs(dx(it)) ~= 0,
@ -86,10 +88,10 @@ while i<n,
% ic=1; % ic=1;
% end % end
end end
if abs(dx(it))>(2*htol), % if abs(dx(it))>(2*htol),
h1(i)= htol/abs(dx(it))*h1(i); % h1(i)= htol/abs(dx(it))*h1(i);
xh1(i)=x(i)+h1(i); % xh1(i)=x(i)+h1(i);
end % end
try try
fx = feval(func,xh1,varargin{:}); fx = feval(func,xh1,varargin{:});
catch catch

52
matlab/mr_hessian.m
View File

@ -50,7 +50,8 @@ if init,
htol = 1.e-4; htol = 1.e-4;
%h1=max(abs(x),gstep_*ones(n,1))*eps^(1/3); %h1=max(abs(x),gstep_*ones(n,1))*eps^(1/3);
%h1=max(abs(x),sqrt(gstep_)*ones(n,1))*eps^(1/6); %h1=max(abs(x),sqrt(gstep_)*ones(n,1))*eps^(1/6);
h1=max(abs(x),sqrt(gstep_)*ones(n,1))*eps^(1/4); % h1=max(abs(x),sqrt(gstep_)*ones(n,1))*eps^(1/4);
h1=options_.gradient_epsilon*ones(n,1);
return, return,
end end
func = str2func(func); func = str2func(func);
@ -102,7 +103,8 @@ while i<n,
icount = 0; icount = 0;
h0=h1(i); h0=h1(i);
while (abs(dx(it))<0.5*htol | abs(dx(it))>(2*htol)) & icount<10 & ic==0, % while (abs(dx(it))<0.5*htol | abs(dx(it))>(2*htol)) & icount<10 & ic==0,
while (abs(dx(it))<0.5*htol) & icount<10 & ic==0,
%while abs(dx(it))<0.5*htol & icount< 10 & ic==0, %while abs(dx(it))<0.5*htol & icount< 10 & ic==0,
icount=icount+1; icount=icount+1;
%if abs(dx(it)) ~= 0, %if abs(dx(it)) ~= 0,
@ -127,21 +129,21 @@ while i<n,
fx=1.e8; fx=1.e8;
end end
end end
if abs(dx(it))>(2*htol), % if abs(dx(it))>(2*htol),
h1(i)= htol/abs(dx(it))*h1(i); % h1(i)= htol/abs(dx(it))*h1(i);
xh1(i)=x(i)+h1(i); % xh1(i)=x(i)+h1(i);
try % try
[fx, ffx]=feval(func,xh1,varargin{:}); % [fx, ffx]=feval(func,xh1,varargin{:});
catch % catch
fx=1.e8; % fx=1.e8;
end % end
while (fx-f0)==0, % while (fx-f0)==0,
h1(i)= h1(i)*2; % h1(i)= h1(i)*2;
xh1(i)=x(i)+h1(i); % xh1(i)=x(i)+h1(i);
[fx, ffx]=feval(func,xh1,varargin{:}); % [fx, ffx]=feval(func,xh1,varargin{:});
ic=1; % ic=1;
end % end
end % end
it=it+1; it=it+1;
dx(it)=(fx-f0); dx(it)=(fx-f0);
h0(it)=h1(i); h0(it)=h1(i);
@ -171,7 +173,7 @@ while i<n,
else else
ff1=ffx; ff1=ffx;
end end
if hflag, % two point based derivatives % if hflag, % two point based derivatives
xh1(i)=x(i)-h1(i); xh1(i)=x(i)-h1(i);
% c=mr_nlincon(xh1,varargin{:}); % c=mr_nlincon(xh1,varargin{:});
% ic=0; % ic=0;
@ -193,9 +195,9 @@ while i<n,
% [f1(:,i), ff1]=feval(func,xh1,varargin{:}); % [f1(:,i), ff1]=feval(func,xh1,varargin{:});
% end % end
ggh(:,i)=(ff1-ff_1)./(2.*h1(i)); ggh(:,i)=(ff1-ff_1)./(2.*h1(i));
else % else
ggh(:,i)=(ff1-ff0)./h1(i); % ggh(:,i)=(ff1-ff0)./h1(i);
end % end
xh1(i)=x(i); xh1(i)=x(i);
if hcheck & htol<1, if hcheck & htol<1,
htol=min(1,max(min(abs(dx))*2,htol*10)); htol=min(1,max(min(abs(dx))*2,htol*10));
@ -209,11 +211,11 @@ h_1=h1;
xh1=x; xh1=x;
xh_1=xh1; xh_1=xh1;
if hflag, % if hflag,
gg=(f1'-f_1')./(2.*h1); gg=(f1'-f_1')./(2.*h1);
else % else
gg=(f1'-f0)./h1; % gg=(f1'-f0)./h1;
end % end
if hflag==2, if hflag==2,
gg=(f1'-f_1')./(2.*h1); gg=(f1'-f_1')./(2.*h1);

15
matlab/newrat.m
View File

@ -82,6 +82,7 @@ else
hhg=hh; hhg=hh;
igg=inv(hh); igg=inv(hh);
end end
H = igg;
disp(['Gradient norm ',num2str(norm(gg))]) disp(['Gradient norm ',num2str(norm(gg))])
ee=eig(hh); ee=eig(hh);
disp(['Minimum Hessian eigenvalue ',num2str(min(ee))]) disp(['Minimum Hessian eigenvalue ',num2str(min(ee))])
@ -105,7 +106,7 @@ while norm(gg)>gtol & check==0 & jit<nit,
bayestopt_.penalty = fval0(icount); bayestopt_.penalty = fval0(icount);
disp([' ']) disp([' '])
disp(['Iteration ',num2str(icount)]) disp(['Iteration ',num2str(icount)])
[fval x0 fc retcode] = csminit(func0,xparam1,fval0(icount),gg,0,igg,varargin{:}); [fval x0 fc retcode] = csminit(func0,xparam1,fval0(icount),gg,0,H,varargin{:});
if igrad, if igrad,
[fval1 x01 fc retcode1] = csminit(func0,x0,fval,gg,0,inx,varargin{:}); [fval1 x01 fc retcode1] = csminit(func0,x0,fval,gg,0,inx,varargin{:});
if (fval-fval1)>1, %(fval0(icount)-fval), if (fval-fval1)>1, %(fval0(icount)-fval),
@ -116,7 +117,7 @@ while norm(gg)>gtol & check==0 & jit<nit,
fval=fval1; fval=fval1;
x0=x01; x0=x01;
end end
if (fval0(icount)-fval)<1.e-2*(gg'*(igg*gg))/2 & igibbs, if (fval0(icount)-fval)<1.e-2*(gg'*(H*gg))/2 & igibbs,
if length(find(ig))<nx, if length(find(ig))<nx,
ggx=ggx*0; ggx=ggx*0;
ggx(find(ig))=gg(find(ig)); ggx(find(ig))=gg(find(ig));
@ -126,6 +127,9 @@ while norm(gg)>gtol & check==0 & jit<nit,
[fvala x0 fc retcode] = csminit(func0,x0,fval,ggx,0,iggx,varargin{:}); [fvala x0 fc retcode] = csminit(func0,x0,fval,ggx,0,iggx,varargin{:});
end end
[fvala, x0, ig] = mr_gstep(0,x0,func0,htol,varargin{:}); [fvala, x0, ig] = mr_gstep(0,x0,func0,htol,varargin{:});
% if length(find(ig))==0,
% [fvala, x0, ig] = mr_gstep(0,x0,func0,htol/10,varargin{:});
% end
nig=[nig ig]; nig=[nig ig];
if (fval-fvala)<gibbstol*(fval0(icount)-fval), if (fval-fvala)<gibbstol*(fval0(icount)-fval),
igibbs=0; igibbs=0;
@ -160,6 +164,8 @@ while norm(gg)>gtol & check==0 & jit<nit,
x(:,icount+1)=xparam1; x(:,icount+1)=xparam1;
fval0(icount+1)=fval; fval0(icount+1)=fval;
%if (fval0(icount)-fval)<ftol*ftol & flagg==1;, %if (fval0(icount)-fval)<ftol*ftol & flagg==1;,
mr_gstep(1,x);
mr_hessian(1,x);
if (fval0(icount)-fval)<ftol, if (fval0(icount)-fval)<ftol,
disp('No further improvement is possible!') disp('No further improvement is possible!')
check=1; check=1;
@ -228,6 +234,7 @@ while norm(gg)>gtol & check==0 & jit<nit,
end end
end end
end end
disp(['Gradient norm ',num2str(norm(gg))]) disp(['Gradient norm ',num2str(norm(gg))])
ee=eig(hh); ee=eig(hh);
disp(['Minimum Hessian eigenvalue ',num2str(min(ee))]) disp(['Minimum Hessian eigenvalue ',num2str(min(ee))])
@ -237,7 +244,9 @@ while norm(gg)>gtol & check==0 & jit<nit,
disp(['Elapsed time for iteration ',num2str(t),' s.']) disp(['Elapsed time for iteration ',num2str(t),' s.'])
g(:,icount+1)=gg; g(:,icount+1)=gg;
save m1.mat x hh g hhg igg fval0 nig % H = bfgsi(H,g(:,end)-g(:,end-1),x(:,end)-x(:,end-1));
H = igg;
save m1.mat x hh g hhg igg fval0 nig H
end end
end end