% Units
% L denotes liter
% M denotes molar which is mole per liter
% nM denotes nano molar

x0 = [0.01; .1; 0; 0];   % nM
D = [ -1   1   1    
      -1   1   0    
       1  -1  -1   
       0   0   1 ];
k = [10 2 0.02];  % [{nM.sec}^{-1} sec^{-1} sec^{-1}]

V = 1e-12;                     % Volume of the compartment  (L)
NA = 6.02214199e23;       % Avogadro's number (molecules per mole)

ssz = NA*V;                    % system size (molecules per M)
ssz = ssz*1e-9;                % system size (molecules per nM)
              
n0 = x0*ssz;                   % molecules

L = -D.*(D < 0);               % Stoich coeffs of the reactants
K = sum(L);                     % molecularity
c = k.*prod(factorial(L))./ssz.^(K-1);
                                     % stochastic rate constant 
tf = 5000;                      % simulation stop time in seconds 
runs = 2;                       % number of realizations(simulation runs)