echo off;
hold off;
close;
c=3;
while c~=0;
 clc, home;
 disp('                       Welcome to the equi-angular spiral by Woody Burley.');
 disp(char(13))
 disp('you may do the following:')
 disp('   0) quit')
 disp('   1) view the deffinition and parametric equations of the equiangular spiral')
 disp('   2) illustrate the definition by choosing an alpha and t to graph tangents and radial')
 disp('      vectors')
 disp('   3) view some real world occurances of the equiangular spiral and some of the history')
 disp('      behind it')
 disp(char(13))
 disp('whenever a graph is showing, just press any key to proceed.')
 c=input('type the # of the item you would like to do >>');
 
 if c==1
  close;
  clc, home;
  disp('         Definition of the equiangular spiral')
  disp(char(13))
  disp('Definition: from any point P on a spiral, draw a line L toward the center of the spiral')
  disp('(usually the origin).  The angle formed by L and the tangent at P is alpha.  If alpha is')
  disp('constant for all points P, the curve is equiangular.')
  disp('the parametric equation for an equiangular spiral in relation to time is:')
  disp('  x=e^(t*cot(alpha))*cos(t),')
  disp('  y=e^(t*cot(alpha))*sin(t).')
  disp('the cartesian equation is:')
  disp('  x^2+y^2=e^(arctan(y/x)*cot(alpha))')
  disp(char(13))
  disp('Some properties of this curve:')
  disp('  when alpha > pi/2, the curve spirals outward towards infinity...')
  disp('  when alpha < pi/2, the curve spirals towards the origin...')
  disp('  when alpha=pi/2, the curve is a circle.')
  disp(char(13));
  disp('  if you draw a radial line L towards infinity, the length from one intersection of L and the')
  disp('  curve, and the length to the next intersection is a geometric progression.... the ratio of')
  disp('  the lengths adjacent to each other is a constant value dependant on alpha.')
  disp('  a few graphs follow; in all of them, t ranges from 0 to 2*pi, the radial vector is green,')
  disp('  and the tangent is red.')
  disp(char(13));
  disp('press any key when you are ready to view the graphs...')
  pause;
  e=2.718281828;
  t=linspace(0,2*pi,100);
  
  T=19*pi/10;
  alpha=pi/6;
  x=e.^(t.*cot(alpha)).*cos(t);
  y=e.^(t.*cot(alpha)).*sin(t);
  plot(x,y);
  axis equal;
  v=axis;
  hold on;
  xlabel('x-axis');
  ylabel('y-axis');
  x1=e^(T*cot(alpha))*cos(T);
  y1=e^(T*cot(alpha))*sin(T);
  plot(x1,y1,'o');
  hold on;
  x2=e^((T+.1)*cot(alpha))*cos(T+.1);
  y2=e^((T+.1)*cot(alpha))*sin(T+.1);
  x4=e^((T-.1)*cot(alpha))*cos(T-.1);
  y4=e^((T-.1)*cot(alpha))*sin(T-.1);
  slope=(y2-y4)/(x2-x4);
  x3=linspace(v(1),v(2));
  y3=y1+slope*(x3-x1);
  plot(x3,y3,'r');
  hold on;
  slope2=y1/x1;
  xo=linspace(v(1),v(2));
  yo=y1+slope2*(xo-x1);
  plot(xo,yo,'g');
  axis(v);
  title('x=e^(^t^*^c^o^t^(^a^l^p^h^a^)^)*cos(t), y=e^(^t^*^c^o^t^(^a^l^p^h^a^)^)*sin(t), alpha=pi/6');
  grid on;
  text(x1,y1,'      alpha=pi/6');
  shg;
  pause;
  close;
  
  alpha=pi/4;
  T=7*pi/4;
  x=e.^(t.*cot(alpha)).*cos(t);
  y=e.^(t.*cot(alpha)).*sin(t);
  plot(x,y);
  axis equal;
  v=axis;
  hold on;
  xlabel('x-axis');
  ylabel('y-axis');
  x1=e^(T*cot(alpha))*cos(T);
  y1=e^(T*cot(alpha))*sin(T);
  plot(x1,y1,'o');
  hold on;
  x2=e^((T+.1)*cot(alpha))*cos(T+.1);
  y2=e^((T+.1)*cot(alpha))*sin(T+.1);
  x4=e^((T-.1)*cot(alpha))*cos(T-.1);
  y4=e^((T-.1)*cot(alpha))*sin(T-.1);
  slope=(y2-y4)/(x2-x4);
  x3=linspace(v(1),v(2));
  y3=y1+slope*(x3-x1);
  plot(x3,y3,'r');
  hold on;
  slope2=y1/x1;
  xo=linspace(v(1),v(2));
  yo=y1+slope2*(xo-x1);
  plot(xo,yo,'g');
  axis(v);
  title('x=e^(^t^*^c^o^t^(^a^l^p^h^a^)^)*cos(t), y=e^(^t^*^c^o^t^(^a^l^p^h^a^)^)*sin(t), alpha=pi/4');
  text(x1,y1-20,'       alpha=pi/4')
  grid on;
  shg;
  pause;
  close;
  
  alpha=pi/3;
  T=15*pi/8;
  x=e.^(t.*cot(alpha)).*cos(t);
  y=e.^(t.*cot(alpha)).*sin(t);
  plot(x,y);
  axis equal;
  hold on;
  v=axis;
  xlabel('x-axis');
  ylabel('y-axis');
  x1=e^(T*cot(alpha))*cos(T);
  y1=e^(T*cot(alpha))*sin(T);
  plot(x1,y1,'o');
  hold on;
  x2=e^((T+.1)*cot(alpha))*cos(T+.1);
  y2=e^((T+.1)*cot(alpha))*sin(T+.1);
  x4=e^((T-.1)*cot(alpha))*cos(T-.1);
  y4=e^((T-.1)*cot(alpha))*sin(T-.1);
  slope=(y2-y4)/(x2-x4);
  x3=linspace(v(1),v(2));
  y3=y1+slope*(x3-x1);
  plot(x3,y3,'r');
  hold on;
  slope2=y1/x1;
  xo=linspace(v(1),v(2));
  yo=y1+slope2*(xo-x1);
  plot(xo,yo,'g');
  axis(v);
  title('x=e^(^t^*^c^o^t^(^a^l^p^h^a^)^)*cos(t), y=e^(^t^*^c^o^t^(^a^l^p^h^a^)^)*sin(t), alpha=pi/3');
  text(x1,y1,'    alpha=pi/3')
  grid on;
  shg;
  pause;
  close;
  
  alpha=10*pi/21;
  x=e.^(t.*cot(alpha)).*cos(t);
  y=e.^(t.*cot(alpha)).*sin(t);
  plot(x,y);
  axis equal;
  hold on;
  v=axis;
  x1=e^(T*cot(alpha))*cos(T);
  y1=e^(T*cot(alpha))*sin(T);
  plot(x1,y1,'o');
  hold on;
  x2=e^((T+.1)*cot(alpha))*cos(T+.1);
  y2=e^((T+.1)*cot(alpha))*sin(T+.1);
  x4=e^((T-.1)*cot(alpha))*cos(T-.1);
  y4=e^((T-.1)*cot(alpha))*sin(T-.1);
  slope=(y2-y1)/(x2-x1);
  x3=linspace(v(1),v(2));
  y3=y1+slope*(x3-x1);
  plot(x3,y3,'r');
  hold on;
  slope2=y1/x1;
  xo=linspace(v(1),v(2));
  yo=y1+slope2*(xo-x1);
  plot(xo,yo,'g');
  axis(v);
  xlabel('x-axis');
  ylabel('y-axis');
  title('x=e^(^t^*^c^o^t^(^a^l^p^h^a^)^)*cos(t), y=e^(^t^*^c^o^t^(^a^l^p^h^a^)^)*sin(t), alpha=10*pi/21');
  text(x1-.5,y1+.35,'alpha=10*pi/21')
  grid on;
  shg;
  pause;
  close;
  
  T=5*pi/3;
  alpha=pi/2;
  x=e.^(t.*cot(alpha)).*cos(t);
  y=e.^(t.*cot(alpha)).*sin(t);
  plot(x,y);
  axis equal;
  hold on;
  v=axis;
  xlabel('x-axis');
  ylabel('y-axis');
  x1=e^(T*cot(alpha))*cos(T);
  y1=e^(T*cot(alpha))*sin(T);
  plot(x1,y1,'o');
  hold on;
  x2=e^((T+.1)*cot(alpha))*cos(T+.1);
  y2=e^((T+.1)*cot(alpha))*sin(T+.1);
  x4=e^((T-.1)*cot(alpha))*cos(T-.1);
  y4=e^((T-.1)*cot(alpha))*sin(T-.1);
  slope=(y2-y4)/(x2-x4);
  x3=linspace(v(1),v(2));
  y3=y1+slope*(x3-x1);
  plot(x3,y3,'r');
  hold on;
  slope2=y1/x1;
  xo=linspace(v(1),v(2));
  yo=y1+slope2*(xo-x1);
  plot(xo,yo,'g');
  axis(v);
  title('x=e^(^t^*^c^o^t^(^a^l^p^h^a^)^)*cos(t), y=e^(^t^*^c^o^t^(^a^l^p^h^a^)^)*sin(t), alpha=pi/2');
  text(x1-.5,y1,'alpha=pi/2')
  grid on;
  shg;
  pause;
  close;
  clear x* y* T v t e;
     
 elseif c==2
   clc, home;    
   alpha=input('choose an angle alpha to generate the curve >>');
   cont=1;
  while cont==1
   clc, home;  
   T=input('choose a time T between 0 and 4pi to find a point on the curve >>');
   e=2.718281828;
   t=linspace(0,4*pi,100);
   x=e.^(t.*cot(alpha)).*cos(t);
   y=e.^(t.*cot(alpha)).*sin(t);
   plot(x,y), xlabel('x-axis'), ylabel('y-axis'), title('user generated graph, with radial vector and tangent line'),axis equal;
   v=axis;
   hold on;
   x1=e^(T*cot(alpha))*cos(T);
   y1=e^(T*cot(alpha))*sin(T);
   plot(x1,y1,'o');
   hold on;
   x2=e^((T+.1)*cot(alpha))*cos(T+.1);
   y2=e^((T+.1)*cot(alpha))*sin(T+.1);
   x4=e^((T-.1)*cot(alpha))*cos(T-.1);
   y4=e^((T-.1)*cot(alpha))*sin(T-.1);
   slope=(y2-y4)/(x2-x4);
   x3=linspace(v(1),v(2));
   y3=y1+slope*(x3-x1);
   plot(x3,y3,'r');
   hold on;
   slope2=y1/x1;
   xo=linspace(v(1),v(2));
   yo=y1+slope2*(xo-x1);
   plot(xo,yo,'g');
   grid on,axis(v);
   format compact;
   disp('the point you chose is:')
   x1,y1
   disp('slope of the tangent line (red):')
   slope
   disp('slope of radial vector (green):')
   slope2
   disp('as you can see, the angle between the radial vector and the tangent is')
   disp(alpha)
   disp(char(13))
   disp('press any key')
   pause;
   close;
   cont=input('enter 1 if you want to choose another point on the spiral, or choose 0 to quit >>');
   clear x y x1 y1 x2 y2 x3 y3 x4 y4 xo yo v t T e;
  end 
   clear alpha cont;
   close;
   
 elseif c==3
  clc, home;
  disp('  The equiangular spiral was discovered by the philosopher/mathemetician Rene Descartes')
  disp('in 1638 and was studied by several other mathemeticians since then (James Bernoulli')
  disp('was so fascinated with the spiral, that he had it engraved on his tomb).  Bernoulli')
  disp('was intrigued by its fractal qualities--if you enlarge or shrink it, it remains the')
  disp('same shape.')
  disp('  Another quality of the spiral is that it is the trace of a "pursuit curve." Or, if')
  disp('there were 4 objects placed an equal distance apart on a circle, and each was attracted')
  disp('to the object immediately counterclockwise from it and they all moved at the same speed,')
  disp('the path they would all take is a equiangular spiral.')
  disp('  In the graphic that follows, the lines are the direction that the point sharing their')
  disp('color will take in the next iteration, there is quite a bit of error due to the size of')
  disp('dt (change in time).')
  disp('press any key.')
  pause;
  t=linspace(0, 2*pi, 20);
  e=2.718281828;
  x=e.^(t.*cot(pi/3)).*cos(t);
  y=e.^(t.*cot(pi/3)).*sin(t);
  for i=15:-1:1
   plot(x(1,i),-y(1,i),'ro')
   hold on;
   plot(y(1,i),x(1,i),'go')
   hold on;
   plot(-x(1,i),y(1,i),'ko')
   hold on;
   plot(-y(1,i),-x(1,i),'bo')
   hold on;
   l1=line([x(1,i) y(1,i)],[-y(1,i) x(1,i)]);
   set(l1, 'color', 'r');
   l2=line([y(1,i) -x(1,i)],[x(1,i) y(1,i)]);
   set(l2,'color','g');
   l3=line([-x(1,i) -y(1,i)],[y(1,i) -x(1,i)]);
   set(l3,'color','k');
   l4=line([-y(1,i) x(1,i)],[-x(1,i) -y(1,i)]);
   set(l4,'color','b');
   if i==15
    axis equal, grid on;
   end   
   pause(1);
end
  clear x y e t i;
  pause;
  close;
  clc, home;
  disp('  Of course, that rarely happens in the real world... equiangular spirals that do happen,')
  disp('in nature are things like snail shells where the growth of an organism is proportional to')
  disp('its size...')
  disp('press any key')
  pause;
  e=2.718281828;
  t=linspace(0,10*pi+3*pi/4,200);
  x=e.^(t.*cot(10*pi/21)).*cos(t);
  y=e.^(t.*cot(10*pi/21)).*sin(t);
  T=8*pi+3*pi/4;
  px(1)=e^(T*cot(10*pi/21))*cos(T);
  py(1)=e^(T*cot(10*pi/21))*sin(T);
  T=10*pi+3*pi/4;
  px(2)=e^(T*cot(10*pi/21))*cos(T);
  py(2)=e^(T*cot(10*pi/21))*sin(T);
  plot(x,y,px,py,'b');
  xlabel('x-axis'), ylabel('y-axis'), title('snail... more or less');
  axis equal, v=axis;, axis([v(1)-1 v(2)+1 v(3)-1 v(4)+1]);
  pause;
  clear x y e t T px py;
  close;
 end
end
clear all
t=linspace(0,10*pi,200);
e=2.718281828;
x=e.^(t.*cot(20*pi/41)).*cos(t);
y=3.^(t.*cot(20*pi/41)).*sin(t);
title('good-bye');
comet(x,y);
title('good-bye');
pause(2);
close, clear all, clc, home;