%   this is a demo for some quadric shapes
% 
 
o = [2,2,0];   % origin will be shifted otherwise we won't see the arrows inside our ellipsoid

 arrow3([0,0,0], [1,0,0],'yellow');
 arrow3([0,0,0], [0,1,0],'yellow');
 arrow3([0,0,0], [0,0,1],'magenta');    % this plots the i,j,k, vectors
 
 xmin = -2;
 xmax = 2;       % plotting window
 
 
   %%  Example 1. Elliptic paraboloid
   %
   %   its equation is  z = x^2 + y^2 
   
 [X,Y] = meshgrid(xmin:0.2:xmax); 
   Z = X.^2 + Y.^2;   
          surf(X,Y,Z);        % surface plotting
 
    
     % traces when intersected with x = k
    
   zmin = -4;
   zmax= 4;
   [Z,Y] = meshgrid(zmin:0.2:zmax, xmin:0.2:xmax); 
   
   k = 1;
   
   X = 0*Y + k; 
     surf(X,Y,Z, 'FaceColor',[0.8 1 0.8]);
             
 
   %%  Example 2. "Saddle" that is hyperbolic paraboloid 
   %
   %   its equation is  z = x^2 - y^2 
   
   
   [X,Y] = meshgrid(xmin:0.2:xmax); 
   Z = X.^2 - Y^2;   
 
         surf(X,Y,Z);        % surface plotting
         
           %% another plot, involving contours 
     
     
      ezmeshc('x^2+ 4*y^2',[-2,2,-2,2])
                    
 
      %% Example 3.  Hyperboloid of two sheets 
      %    Z^2 = 
      
      % we will just plot one, the other one is for -Z 
      
     % Z = 1 + (X^.2 + Y^.2); 
        
      ezmeshc('sqrt(1+x^2 +y^2)',[-2,2,-2,2])
         hold on
      ezmeshc('-sqrt(1+x^2 +y^2)',[-2,2,-2,2])   
          hold off
      
        % you will get a double cone if  you drop 1+  from above 
      
      
      
      
      
      
      
      
      
      
      
      
 
      
      
      %% Example 3: Ellipsoid with uneven axes 
      %   probably hard! 
      
 % 
 %     its equation is 
 %     5x^2 - 4xy + 5y^2 + z^2/4 = 4 
 %
 
 
 [X,Y] = meshgrid(xmin:0.2:xmax); 
   Zero = 0*X;
 Z = 4*sqrt( max(1 - (X-Y).^2 - (X+Y).^2/4, Zero) );   % have to avoid taking sqrt(negative)
 
         surf(X,Y,Z);        % surface plotting
         
        
 % this is a nicer plot using spherical coordinates
 
 d = pi/48;
         
 [phi, theta] = meshgrid(-pi/2:d:pi/2,  (0:d:2*pi));
 
   x1 =  cos(phi).* cos(theta);
   y1 =  cos(phi).* sin(theta);
   z1 =  sin(phi); 
   
   % now our surface is actually transformed from x1,y1,z1 
   
   x = (x1 + 2*y1)/2;
   y = (-x1 + 2*y1)/2;
   z = 4*z1;
 
   surf(x,y,z);