function [ sspoint ] = GeneSwitch( startpoint )
% uses crude ODE solver to identify steady-state concentrations of simple,
% mutually repressive genetic switch
% takes a vector (x,y) of initial state
% returns a vector of steady state

residual = 1;
% simulation and model parameters
dt = 0.1;
mu = 20;

x = startpoint(1); y = startpoint(2);
% loop until converged
while abs(residual) > 0.00001
    % update concentrations
    nx = x + dt*(mu/(1+y^2) - x);
    ny = y + dt*(mu/(1+x^2) - y);
    residual = nx-x;
    x = nx; y = ny;
end
sspoint = [x y];

end

