##  M. Turhan Coban
# f(x()) multifunction multi x
function df = dfi (f,x,dx)
if nargin<3 
    dx=0.001;
end
fi=f(x);
n=length(fi);
m=length(x);
  for i=1:n
    for j=1:m
      df(i,j)=0.0;
    end
  end
  for i=1:n
    for j=1:m
      xi=x;
      xi(j)=xi(j)-4.0*dx;
      ff=f(xi);
      df(i,j)=df(i,j)+3.0*ff(i);
      xi=x;
      xi(j)=xi(j)-3.0*dx;
      ff=f(xi);
      df(i,j)=df(i,j)-32.0*ff(i);
      xi=x;
      xi(j)=xi(j)-2.0*dx;
      ff=f(xi);
      df(i,j)=df(i,j)+168.0*ff(i);
      xi=x;
      xi(j)=xi(j)-dx;
      ff=f(xi);
      df(i,j)=df(i,j)-672.0*ff(i);
      xi=x;
      xi(j)=xi(j)+dx;
      ff=f(xi);
      df(i,j)=df(i,j)+672.0*ff(i);
      xi=x;
      xi(j)=xi(j)+2.0*dx;
      ff=f(xi);
      df(i,j)=df(i,j)-168.0*ff(i);
      xi=x;
      xi(j)=xi(j)+3.0*dx;
      ff=f(xi);
      df(i,j)=df(i,j)+32.0*ff(i);
      xi=x;
      xi(j)=xi(j)+4.0*dx;
      ff=f(xi);
      df(i,j)=df(i,j)-3.0*ff(i);
      df(i,j)=df(i,j)/(840.0*dx);
    end
  end
endfunction
