 

clear
close

disp('     THE GRADIENT VECTOR FIELD    ')

disp('            f(x,y) =  xy - x^3/3  ')
disp('         gradf = ( y - x^2, x )  ')

 f= inline(' x.*y - (x.^3)/3', 'x', 'y');
 fx = inline(' y - x.^2', 'x', 'y');
 fy = inline(' x', 'x', 'y');
 
  x = linspace( -2,2, 101);
  y = linspace( -1.5, 1.5, 81); 

  [X,Y] = meshgrid(x,y);

  surf(X,Y,f(X,Y)); shading interp;
  hold on
  levels = linspace(-6,6, 25);
  contour3(X,Y,f(X,Y), levels,'k')
  view(10,50)
  rotate3d on
  pause



 figure
 contour(X,Y,f(X,Y), levels)
 hold on

  xx = linspace(-2,2,21);
  yy = linspace( -1.5, 1.5, 16);
  [XX,YY] = meshgrid(xx,yy);
  U = fx(XX,YY); V = fy(XX,YY);
  quiver(XX,YY, U,V);
  hold off

  pause
  close all

disp('  g(x,y) =  2*exp(-(x-1)^2 -(y-1)^2) -exp(-2*(x+.3)^2 -3*(x+.5)^2) ')

clear U V  levels

g1 = '2*exp(-(x-1).^2 - (y-1).^2)';
g2 = ' -exp(-2*(x+.3).^2 - 3*(y+.5).^2)';
g = inline([g1,g2], 'x', 'y');

surf(X,Y,g(X,Y)); shading interp;
hold on
levels = linspace(-1,2,21);
contour3(X,Y,g(X,Y), levels, 'k')
rotate3d on
pause

figure

contour(X,Y,g(X,Y), levels)
hold on


 g1x = ' -4*(x-1).*exp(-(x-1).^2 - (y-1).^2)';
 g2x = '+4*(x+.3).*exp(-(x+.3).^2 -3*(y+.5).^2)';
 gx = inline([g1x,g2x], 'x', 'y');

 g1y = '-4*(y-1).*exp(-(x-1).^2 -(y-1).^2)';
 g2y = '+6*(y+.5).*exp(-2*(x+.3).^2 -3*(y+.5).^2)';
 gy = inline([g1y,g2y], 'x', 'y');

  U = gx(XX,YY);
  V = gy(XX,YY);
  quiver(XX,YY,U,V)
  axis equal
  pause
  close all
