dynprog.m

dynamic programming 源代码

function [p_opt, fval]=dynprog(x, DecisFun, ObjFun, TransFun) 
% Input x column state variable vector 
% DecisFun(k, x) decision variable at state x of stage k 
% ObjFun(k, x,u) object function at state x and decision u of stage k 
% TranFun(k,n,u) transition function 
% Output p_opt=[order; optimal route, optimal Decision; function value] 
% fval optimal function value 
k = length(x(1,:); 
f_opt = nan*ones(size(x)); 
d_opt = f_opt; 
t_vubm = inf*ones(size(x)); 
x_isnan = ~isnan(x); 
t_vub=inf; 
tmp1 = find(x_isnan(:,k)); 
tmp2 = length(tmp1); 
for i = 1:tmp2 
u = feval(DecisFun, k, x(tmp1(i),k)); 
tmp3 =length(u); 
for j=1:tmp3 
tmp=feval(ObjFun,k,x(tmp1(i),k),u(j)); 
if tmp<=t_vub, 
f_opt(i,k)=tmp; 
d_opt(i,k)=u(j); 
t_vub=tmp; 
end 
end 
end 
for ii=k-1:-1:1 
tmp10=find(x_isnan(:,ii)); 
tmp20=length(tmp10); 
for i=1:tmp20 
u=feval(DecisFun, ii, x(tmp10(i),ii)); 
tmp30 =length(u); 
for j=1:tmp30 
tmp00=feval(ObjFun, ii, x(tmp10(i),ii),u(j)); 
tmp40=feval(TransFun, ii,x(tmp10(i),ii),u(j)); 
tmp50=x(:,ii+1)-tmp40; 
tmp60=find(tmp50==0); 
if ~isempty(tmp60) 
tmp00=tmp00+f_opt(tmp60(1),ii+1); 
if tmp00<=t_vubm(i,ii) 
f_opt(i,ii)=tmp00; 
d_opt(i,ii)=u(j); 
t_vubm(i,ii)=tmp00; 
end 
end 
end 
end 
end 
p_opt=[]; 
tmpx=[]; 
tmpd=[]; 
tmpf=[]; 
tmp0=find(x_isnan(:,1)); 
tmp01=length(tmp0); 
fval=f_opt(tmp01,1); 
for i=1:tmp01 
tmpd(i)=d_opt(tmp0(i),1); 
tmpx(i)=x(tmp0(i),1); 
tmpf(i)=feval(ObjFun,1,tmpx(i),tmpd(i)); 
p_opt(k*(i-1)+1, [1,2,3,4])=[1,tmpx(i),tmpd(i),tmpf(i)]; 
for ii =2:k 
tmpx(i)=feval(TransFun, ii-1,tmpx(i),tmpd(i)); 
tmp1=x(:,ii)-tmpx(i); 
tmp2=find(tmp1==0); 
if ~isempty(tmp2) 
tmpd(i)=d_opt(tmp2(1),ii); 
end; 
tmpf(i)=feval(ObjFun,ii,tmpx(i),tmpd(i)); 
p_opt(k*(i-1)+ii,[1,2,3,4])=[ii,tmpx(i),tmpd(i),tmpf(i)]; 
end 
end