clear
clc
p=1;
T0=125;
s_0=100;
P=[1 0 1 0 0 1 1 0 1 0 0 1;0 1 1 0 0 1 0 1 1 0 0 1;0 0 1 1 0 1 0 0 1 1 0 1;0 0 1 0 1 1 0 0 1 0 1 1];
TT=1;
num_T=10;%input('num_T');
while(TT<num_T)
s0=s_0;
q=poissrnd(9.000,4,12)/48;
u=18.2222/48;
x=rand(1,4);count=0;count_p=0;
while (sum(x/sum(x)>2/T0)<4)
x=rand(1,4);
end
x=x/sum(x);
x0=x;
while (count<35&count_p<p)
x=4*x.*(1-x);
x=x/sum(x);
tt=round(x*T0);
t(TT,:)=tt;
t_c=round(x0*T0);%t(1,:);
pijk=0;qijk=0;
for i=1:4
for j=1:4
for k=1:3
pijk=ceil(P(i,3*(j-1)+k)*u*tt(i)+pijk);
qijk=ceil(q(i,3*(j-1)+k)*tt(i)+qijk);
pijk_c=ceil(P(i,3*(j-1)+k)*u*t_c(i)+pijk);
qijk_c=ceil(q(i,3*(j-1)+k)*t_c(i)+qijk);
end
end
end
clear tt
s(TT)=qijk-pijk;
s(TT)=(s(TT)+abs(s(TT)))/2;
s_c(TT)=qijk_c-pijk_c;
s_c(TT)=(s_c(TT)+abs(s_c(TT)))/2;
if s0>s(TT)
s0=s(TT);
x0=x;
else
%s(TT)=s0;
count_p=count_p+1;
end
end
count=count+1;
%qijk
%pijk
TT=TT+1;
end
%t
%t_c
s
%s_c