function Jf = numjaco(f,m,x,n)
% f is a function with m components,
% x is a vector with n components,
% the result is an m by n matrix.
Jf = ones(m,n);   h = 10^(-4);
for j =1:n
    e = zeros(1,n); e(j) = 1;
    Jf(:,j) = (feval(f,x+h*e)-feval(f,x-h*e) )/(2*h);
end;
% FILE  numjaco.m  ends.
