function ss_branch();
alpha=[12 8 3 2 0];
beta=[10 3 5 6 0];
g=zeros(5,5);
h=zeros(5,5);
g(1,3)=-0.8;
g(1,5)=1;
g(2,1)=0.5;
g(3,2)=0.75;
g(4,1)=0.5;
h(1,1)=0.5;
h(2,2)=0.75;
h(3,3)=0.5;
h(3,4)=0.2;
h(4,4)=0.8;
p0=[0.2 0.3 4.2 0.01 1];
tspan = [0:0.25:5];
[t,p] = ode15s(@ode_ssys, tspan, p0,[],alpha,beta,g,h);
figure(1);
plot(t,p(:,1),'rs-');; 
hold on;
plot(t,p(:,2),'gd-');
plot(t,p(:,3),'ko-');
plot(t,p(:,4),'mx-');
legend('X1','X2','X3','X4');
xlabel('Time')
xlabel('Concentration')


function dp = ode_ssys(t,X,alpha,beta,g,h);
dp = zeros(5,1);   
for i=1:5   
  vplus=alpha(i); 
  for j=1:5
     vplus=vplus*X(j)^g(i,j);
  end
  vminus=beta(i);  
  for j=1:5
     vminus=vminus*X(j)^h(i,j);
  end
  dp(i)=vplus-vminus;
end
