function [the_out_exact,P_exact,E_exact,tao_the]=MHDOPP12exact(total_node_num,node,B,C,Icood,signs,Effb,th_out,E_avg,equal)
%%%%%%%outge with a specific buffer conbination
%%%will call MHDOPP13.m, very slow
L=length(node)-1;
for l=L-2:-1:0
Bvector(L-1-l)=(B+1)^l;
end
clear state



%%%%%%%%%%%exact solution
if equal==0
for i=1:(L-1)
    if i==1
        for j=1:B+1
            state(j,1)=j-1;
        end
    end
    if i>1
        state2=[];
        size_state=size(state);
        for i1=0:B
            for j=1:size_state(1,1)
                state2((i1)*size_state(1,1)+j,:)=[i1 state(j,:)];
            end
        end
        state=state2;
    end
end
clear state2 size_sate i i1 j l itg 

%%%generate P and E book 
for i=1:2^(length(Icood)-1)-1%no zero, SD always available, when only consider the impact from buffer
       clear availablechannels
       temp=dec2bin(i,length(Icood)-1);
       for j=1:length(Icood)-1
       availablechannels(j+1)=bin2dec(temp(j));
       end
       availablechannels(1)=1;
       [P_temp,E_temp]=MHDOPP13(total_node_num,C,Icood,signs,Effb,th_out,E_avg,availablechannels);
       P_all(i,:,:)=P_temp;
       E_all(i,:,:)=E_temp;
end

sizestate=size(state);
Markov_trans=zeros(sizestate(1,1),sizestate(1,1));
Markov_trans1=zeros(sizestate(1,1),sizestate(1,1));
for state_num=1:sizestate(1,1)
    %%%%%%%%%%%from state_num to state
    i_temp=state_num-1;
    for l=L-1:-1:1
        if l>1
        temp=floor(i_temp/(B+1)^(l-1));
        state_temp(L-l)=temp;
        i_temp=i_temp-temp*(B+1)^(l-1);
        elseif l==1
        state_temp(L-l)=i_temp;
        end
    end
    
    availablechannels=ones(1,length(Icood));
    for i=1:length(state_temp)
        if state_temp(i)==0
           for j=1:length(Icood)
               if Icood(j,1)==i+2%from node 3
                  availablechannels(j)=0;
               end
           end
        elseif state_temp(i)==B
           for j=1:length(Icood)
               if Icood(j,2)==i+2
                  availablechannels(j)=0;
               end
           end 
        end
    end
    P=reshape(P_all(sum(availablechannels(2:end).*(2.^([length(Icood)-2:-1:0]))),:,:),total_node_num,total_node_num);
    
    for i=2:length(Icood)%ignore SD here, deal with it at last
        input_system=0;
        if P(Icood(i,1),Icood(i,2))>0
           state_temp2=state_temp;
           if Icood(i,1)==1&Icood(i,2)~=2
              state_temp2(Icood(i,2)-2)=state_temp2(Icood(i,2)-2)+1;
              input_system=1;
           elseif Icood(i,2)==2&Icood(i,1)~=1
              state_temp2(Icood(i,1)-2)=state_temp2(Icood(i,1)-2)-1;
           elseif Icood(i,2)~=2&Icood(i,1)~=1
              state_temp2(Icood(i,1)-2)=state_temp2(Icood(i,1)-2)-1;
              state_temp2(Icood(i,2)-2)=state_temp2(Icood(i,2)-2)+1;     
           end
           Markov_trans(state_num,sum(Bvector.*state_temp2)+1)=P(Icood(i,1),Icood(i,2));
           if input_system==1
           Markov_trans1(state_num,sum(Bvector.*state_temp2)+1)=P(Icood(i,1),Icood(i,2));
           end
        end
    end
    Markov_trans(state_num,state_num)=1-sum(Markov_trans(state_num,:));
    
end


for i2=1:length(Markov_trans)
    Markov_trans(i2,:)=Markov_trans(i2,:)/sum(Markov_trans(i2,:));
end
all_state=random('unif',0,1,1,sizestate(1,1));
all_state=all_state/sum(all_state);
for i=1:1000
all_state=all_state*Markov_trans;
end
i3=1;
clear all_state_bak
for i=1:20
all_state_bak(i3,:)=all_state;i3=i3+1;
all_state=all_state*Markov_trans;
end
for i=2:20
    if i==2&sum(abs(all_state_bak(1,:)-all_state_bak(2,:)))<length(all_state)*10e-3
        all_state=all_state_bak(1,:);
    elseif sum(abs(all_state_bak(1,:)-all_state_bak(i,:)))<length(all_state)*10e-3
        all_state=sum(all_state_bak(1:i-1,:))/(i-1);
        break
    end
end

%%%%%%%%%the_out
the_out_exact=0;
P_exact=zeros(total_node_num,total_node_num);
E_exact=zeros(total_node_num,total_node_num);
for state_num=1:length(all_state)
     %%%%%%%%%%%from state_num to state
    i_temp=state_num-1;
    for l=L-1:-1:1
        if l>1
        temp=floor(i_temp/(B+1)^(l-1));
        state_temp(L-l)=temp;
        i_temp=i_temp-temp*(B+1)^(l-1);
        elseif l==1
        state_temp(L-l)=i_temp;
        end
    end
  
    availablechannels=ones(1,length(Icood));
    for i=1:length(state_temp)
        if state_temp(i)==0
           for j=1:length(Icood)
               if Icood(j,1)==i+2%from node 3
                  availablechannels(j)=0;
               end
           end
        elseif state_temp(i)==B
           for j=1:length(Icood)
               if Icood(j,2)==i+2
                  availablechannels(j)=0;
               end
           end 
        end
    end
    the_out_temp=1;
    for i=1:length(availablechannels)
        if availablechannels(i)==1
           the_out_temp=the_out_temp*(1-exp(-th_out(Icood(i,1),Icood(i,2))));
        end
    end
    the_out_exact=the_out_exact+all_state(state_num)*the_out_temp;
    P_exact=P_exact+all_state(state_num)*reshape(P_all(sum(availablechannels(2:end).*(2.^([length(Icood)-2:-1:0]))),:,:),total_node_num,total_node_num);
    E_exact=E_exact+all_state(state_num)*reshape(E_all(sum(availablechannels(2:end).*(2.^([length(Icood)-2:-1:0]))),:,:),total_node_num,total_node_num);
end

tao_the=0;%average packet delay
for i=1:length(all_state)
    tao_the=tao_the+(i-1)*all_state(i);
end
temp=0;
if equal==0
for i=1:length(all_state)%for direct transmission
    if i==1
    temp=temp+all_state(i)*(Markov_trans(i,i)-(1-exp(-th_out(1,2)))*(1-exp(-th_out(1,3))));
    elseif i==length(all_state)
    temp=temp+all_state(i)*(Markov_trans(i,i)-(1-exp(-th_out(1,2)))*(1-exp(-th_out(2,3))));
    else
    temp=temp+all_state(i)*(Markov_trans(i,i)-(1-exp(-th_out(1,2)))*(1-exp(-th_out(2,3)))*(1-exp(-th_out(1,3))));
    end
end

for i=1:length(all_state)-1
    temp=temp+all_state(i)*Markov_trans(i,i+1);
end
tao_the=tao_the/temp+1;
elseif equal==1
    tao_the=0;
end



end%end of equal==0

if equal==1
    %%%%%%%%%the output
    the_out_exact=0;
    P_exact=zeros(total_node_num,total_node_num);
    E_exact=zeros(total_node_num,total_node_num);

for i2=0:3^(total_node_num-2)-1
    availablechannels=ones(1,length(Icood));
    temp=dec2base(i2,3,total_node_num-2);
    for j1=1:total_node_num-2
        tempdec(j1)=base2dec(temp(j1),3);
        for j2=1:length(Icood)
            if Icood(j2,1)==j1+2&&tempdec(j1)==0
               availablechannels(j2)=0;
            elseif Icood(j2,2)==j1+2&&tempdec(j1)==2
               availablechannels(j2)=0;
            end
        end
    end
    
    [P_temp,E_temp]=MHDOPP13(total_node_num,C,Icood,signs,Effb,th_out,E_avg,availablechannels);
    P_all(i2+1,:,:)=P_temp;
    E_all(i2+1,:,:)=E_temp;
    
    state_equal(i2+1)=1;
    for j1=1:total_node_num-2
        if tempdec(j1)==1
           state_equal(i2+1)=state_equal(i2+1)*(B-1)/(B+1);
        else
          state_equal(i2+1)=state_equal(i2+1)*(1)/(B+1);
        end
    end
    
    the_out_temp=1;
    for i=1:length(availablechannels)
        if availablechannels(i)==1
           the_out_temp=the_out_temp*(1-exp(-th_out(Icood(i,1),Icood(i,2))));
        end
    end
    the_out_exact=the_out_exact+state_equal(i2+1)*the_out_temp;
    P_exact=P_exact+state_equal(i2+1)*reshape(P_all(i2+1,:,:),total_node_num,total_node_num);
    E_exact=E_exact+state_equal(i2+1)*reshape(E_all(i2+1,:,:),total_node_num,total_node_num);
    tao_the=0;%not calculate the average packet delay
end
        
end%end of equal==1

end