function s = Milgrom_free_function(t,w,flag,b) % s = zeros(2,1); s(1) = w(2); s(2) = -sqrt(b*abs(w(1)))*sign(w(1));