function [I] = simpson13(x,f)
  [nrx,ncx]=size(x)
  [nrf,ncf]=size(f)
  if((nrx<>1)|(nrf<>1))then
    error('x or f,or both,not colum vector');
    abort;
  end;
  if((ncx<>ncf))then
    error('x and f are not of the same lebgth');
    abort;
  end;
  //check the size of the vector x and f is odd
  if(modulo(ncx,2)==0)then
    disp(ncx,"list size =")
    error('list size must be an odd number');
    abort;
  end;
  n = ncx;
  xdiff = mtlb_diff(x);
  h = xdiff(1,1);
  I = f(1,1) + f(1,n);
  for j 2:n-1
    if(modulo(j,2)==0) then
      I = I + 4*f(1,j);
    else
      I = I + 2*f(1,j);
    end;
  end;
  I = (h/3.0)*I
// end of function simpson 1/3