clc;clear all;close all;
q=[.02 .02 .02 .02 .02  .1 .1 .1 .1     .01 .01 .01 .01 .01 .01     .02 .02 .02 .02     .04 .04 .04     .04 .04 .04 .04     .05 .05 .05     .08 .12 .12];
C=[12 12 12 12 12      20 20 20 20      50  50  50  50  50  50      76  76  76  76     100 100 100      155 155 155 155    197 197 197     350 400 400];

 X=ones(1,3406);
 for i=1:size(X,2)
     X(1,i)=i-1;
 end

 COPT=zeros(3406,5,32);
for k=1:size(COPT,3)
    COPT(:,1,k)=X';
end

 COPT(1,2,1)=1;

for k=1:size(COPT,3)
    
    if k==1
        for i=1:size(COPT,1)
            COPT(i,3,k)=COPT(i,1,k)-C(k);
            if COPT(i,3,k)<=0 
                COPT(i,4,k)=1;
            end
            COPT(i,5,k)=(1-q(k))*COPT(i,2,k) + q(k)*COPT(i,4,k);
        end
     
    else
        for i=1:size(COPT,1)
            COPT(i,3,k)=COPT(i,1,k)-C(k);
            COPT(i,2,k)=COPT(i,2,k-1);
            if COPT(i,3,k)<=0 
                COPT(i,4,k)=1;
            else
                COPT(i,4,k)= COPT(i-C(k),5,k-1);
            end

            COPT(i,2,k)=COPT(i,5,k-1);
            COPT(i,5,k)=(1-q(k))*COPT(i,2,k) + q(k)*COPT(i,4,k);
        end
    end
end
fa=COPT(:,5,32);
%a=[0 12 20 24 32 36 40 44 48 50 52 56 60 76 80 100 120 140 160 180 200 250 300 350 400 450 500 550 600 650 700 750 800 850 900 950 1000 1100 1200 1300 1400 1500 1600  1700 1800 1900 2000 2200 2400 2600 2800 3000 3200 3400 3405];
%b=[1 .7636 .7398 .6344 .6334 .6227 .6226 .6051 .6047 .6047 .5904 .5886 .5886 .5792  .5579 .5476];