function [rec, ratio, occ] = GillespieSwain(rateL, rateKR, rateKT, nbind, sugar)
% performs Gillespie simulation on Monod model
% takes L, KR, KT rate parameters; also "expt" (0, 1, or 2: number of sugar binding sites) and sugar copy number
% returns recorded snapshots of the system (rec) and ratio of bound state
% occupancy (ratio)

% initial copy numbers
  x = [1 0 0 0 0 0 sugar];
% simulation parameters
  t = 0; endt = 100;
  N = 7; M = 10;
  a = zeros(M,1);
  rec = zeros(endt,N);
  occ = zeros(1,N);
  timelast = 0;
  
  % stoichiometry and required reactant matrices
  stoichiometry = [ -1 1 0 0 0 0 0;
      1 -1 0 0 0 0 0;
      -1 0 1 0 0 0 -1;
      1 0 -1 0 0 0 1;
      0 0 -1 1 0 0 -1;
      0 0 1 -1 0 0 1;
      0 -1 0 0 1 0 -1;
      0 1 0 0 -1 0 1;
      0 0 0 0 -1 1 -1;
      0 0 0 0 1 -1 1];
  reactants = -stoichiometry.*(stoichiometry < 0);
  % reaction rates
  c = [rateL 1 nbind rateKR nbind-1 2*rateKR nbind rateKT nbind-1 2*rateKT];
  
  while(t < endt)
      % calculation reaction propensities
      h = ones(M,1);
      for j=1:M
          for k=1:N
              h(j) = h(j)*CombinatoricConcentration(x(k), reactants(j,k));
          end
      end
      
      % translate to expected timings
      az0 = 0;
      for j=1:M
          a(j) = h(j)*c(j);
          az0 = az0 + a(j);
      end
      
      % random number draw
      r1 = rand(); r2 = rand();
      tau = (1/az0)*log(1/r1);
      
      % locate reaction corresponding to draw
      tmp = r2*az0; 
      cumsuma = cumsum(a);
      reaction = -1;
      for j=1:M
          if cumsuma(j) > tmp && reaction == -1
              reaction = j;
          end
      end

      % update records and system state
      timenow = floor(t+tau);
      for tmpt=timelast:timenow
        rec(tmpt+1,:) = x;
      end
      occ = occ+tau*x;
      timelast = timenow;
      x = x + stoichiometry(reaction,:);
      t = t + tau;
  end
 
  % calculate binding ratio
  ratio = sum(rec(:,1)+rec(:,3)+rec(:,4))/sum(rec(:,1)+rec(:,3)+rec(:,4)+rec(:,2)+rec(:,5)+rec(:,6));
      
  end