%
%  w05pflowA.m -- (djm:  05 feb 2010)
%

%  initialize workspace
clear;  close all

%  set x & y-vectors &  grid matrices ("help meshgrid")
L = 2.0;  ds = 2*L/101;
x = -L:ds:+L;  y = -L:ds:+L;
[xx,yy] = meshgrid(x,y);  
zz = xx + i*yy;  rr2 = xx.^2 + yy.^2;

%  menu
flo = [];  
while (isempty(flo)==1)
  disp(' ')
  disp('  flow choices:')
  disp('    1) uniform flow')
  disp('    2) line vortex')
  disp('    3) source/sink')
  disp('    4) image sources')
  disp(' ')

  flo = input(' choice:  ');
end

%  complex potential
switch(flo)

%  uniform flow
  case{1}
    L0 = 1;  A = 1.0;  ang = pi/3;
    Phi = A*exp(-i*ang)*zz + 10*i;
    ftitle = ['\bf uniform flow'];
    cI = 40;  cR = 10;
	
%  line vortex
  case{2}
    L0 = 1;  lambda = 1;;
    Phi = -i*(lambda/2/pi)*log(zz);
    ftitle = ['\bf line vortex'];
    cI = 40;  cR = 10;
	
%  source/sink
  case{3}
    L0 = 1;  Q = 1;
    Phi = (Q/2/pi)*log(zz);
    ftitle = ['\bf source/sink'];
    cI = 40;  cR = 10;
	
%  source near wall	
  case{4}
    L0 = 1;  Q = 1;
    Phi = (Q/2/pi)*log(zz.^2-1);
    ftitle = ['\bf image sources'];
    cI = 40;  cR = 10;	
end

%  plot level curves of Phi

clf
contour(x,y,imag(Phi),cI)
caxis([min(imag(Phi(:))) max(imag(Phi(:)))])
colormap(jet);  colorbar
hold on

contour(x,y,real(Phi),cR,'g--')
axis equal;  axis([-1 1 -1 1]*L)

%  label plot

title(ftitle)
xlabel('\bf Re-axis')
ylabel('\bf Im-axis')

text(-0.95*L0*L,0.9*L0*L, ...
    '\bf contours:  streamlines in colour, velocity potential in green')