dynamic.m

清华大学运筹学课件

% the programm is with the dynamic programming use the recurisive method for the last to first
% this is the main function of the method

function[p_opt,fval]=dynprog(x,DecisFun,ObjFun,TransFun)

% the function is to solve the dynamic example in the textbook
% x is the situation variant and its column number represent the stage situation
% subfunction DecisFun(k,x) is to solve the decision variant of k stage variant x 
% subfunction ObjFun(k,x,u) is to stage index function
% subfunction TransFun(k,x,u) is the stage transformation function,u is the corresponding decision variant
% p_opt has four output,the first is the number of the stage,the second is the optimal road of decision
% the third is the optimal stategies of the decision ,the forth is the index function group.
% fval is a column vector,is to represent the optimal value correspend to  the initial stage is x 
% 

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;

%%%%%%%%%%%%%%%%%

% to caculate the teminate values

tmp1=find(x_isnan(:,k));
tmp2=length(tmp1);

for i=1:tmp2
     u=feval(DecisFun,k,x(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
%%%%%%%%%%%%%%%%

% recurisive

for ii=k-1:-1:1
       tmp10=find(x_isnan(:,ii));
       tmp20=length(tmp10);
  for i=1:tmp20
     u=feval(DecisFun,ii,x(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),
         if nargin<5
           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
end

fval=f_opt(tmp1,1);

%%%%%%%%%%%%%%%%%%

% to write the index and parameter or result

p_opt=[];
tmpx=[];
tmpd=[];
tmpf=[];
tmp0=find(x_isnan(:,1));
tmp01=length(tmp0);
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