%Flujo incompresible bidimensional de un sumidero
% externo a un disco impermeable.
% Ocatve: usar " "
% Matlab: usar ' '

%Parametros
%Radio del disco:
a = 1 ;
%Coordenada y magnitud de la fuente:
b = 1.5*a ;
m = 1 ;
%Coordenada del sumidero imagen (magnitud m):
c = a^2/b ;
%(fuente imagen en O, de magnitud m)

Npx = 500 ;
Npy = 500 ;
xmin = -2*a ;
xmax = 2*a ;
ymax = 1.7*a ;
ymin = -ymax ;
levelspsi = -.5:.025:.5 ;
levelsphi = -.5:.025:.5 ;

x = linspace(xmin, xmax, Npx) ;
y = linspace(ymin, ymax, Npy) ;
[xx, yy] = meshgrid(x,y) ;
z = xx + i*yy ;

t = 0:.01:2*pi ;


w = (m/(2*pi))*log( z./((z-b).*(z-c)) ) ;
psi = imag(w) ;
phi = real(w) ;

figure
contour(xx,yy,psi, levelspsi, "g", "linewidth", 1.5)
hold
contour(xx,yy,phi, levelsphi, "b", "linewidth", 1.5)
legend('psi=cte', 'phi=cte')
grid
axis("equal")
set(gca(), "fontsize", 18)

figure
contour(xx,yy,psi, levelspsi, "g", "linewidth", 1.5)
hold
contour(xx,yy,phi, levelsphi, "b", "linewidth", 1.5)
fill(a*cos(t), a*sin(t), [.9 .9 .95])
legend('psi=cte', 'phi=cte')
grid
axis("equal")
set(gca(), "fontsize", 18)