function grad_ascent(x,y,z,px,py,N,mu,xstart,ystart)

xga(1)= xstart;
yga(1)= ystart;
zga(1)=func(xga(1),yga(1));

for i=1:N
  gradx = ( func(xga(i)+eps,yga(i))-func(xga(i),yga(i)) )/eps;
  grady = ( func(xga(i),yga(i)+eps)-func(xga(i),yga(i)) )/eps;
  xga(i+1) = xga(i) + mu*gradx;
  yga(i+1) = yga(i) + mu*grady;
  zga(i+1)=func(xga(i+1),yga(i+1));
end

hold off
contour(x,y,z,10)
hold on
quiver(x,y,px,py)
hold on
plot(xga,yga)
S = sprintf('Gradiant Ascent: N = %d, Step Size = %f',N,mu);
title(S)
xlabel('x axis')
ylabel('yaxis')

