Fluido alrededor de un obstáculo circular (Grupo 20A)
Contenido
1 Introducción
Un fluido es un medio continuo tal que en reposo no transmite tensión tangencial. No es capaz soportar su propio peso, y tiende a adoptar las características del recipiente en el que se encuentra.
Los fluidos pueden ser líquidos o gases, se van a estudiar los fluidos incompresibles, que son aquellos cuyo volumen no disminuye al ejercer una fuerza sobre él, este tipo de fluidos son los líquidos.
En este artículo se estudiará el comportamiento de estos fluidos al fluir alrededor de un obstáculo con forma circular. Por conveniencia se trabajará con coordenadas cilíndricas.
2 Mallado
Se empezará dibujando un mallado que represente los puntos interiores de la región ocupada por un fluido, que será el exterior del círculo unidad. Para llevar a cabo esta gráfica, se tomará un mallado del anillo comprendido entre los radios 1 y 5 y centro el origen.
La finalidad de realizar un mallado es dividir el dominio continuo (la región ocupada por el fluido) en subdominios más pequeños llamados elementos finitos o celdas. Estos elementos forman una malla que cubre toda la región y permite realizar cálculos y simulaciones numéricas en cada uno de ellos.
Para ilustrar que el fluido ocupa el exterior de un circulo, se dibujarán los ejes en el intervalo [math][-5,5]×[-5,5][/math].
rho=linspace(1,5,30); %Se definen rho y theta
th=linspace(0,2*pi,30);
[U,V]=meshgrid(rho,th); %Se crea el mallado
hold on
X=U.*cos(V); %Se parametriza la superficie
Y=U.*sin(V);
Z=0.*U;
mesh(X,Y,0*Z); %Dibujo de la placa
plot(1*cos(th),1*sin(th),'k','lineWidth',1); %Se representa el obstáculo
axis([-5,5,-5,5]); %Se fijan los ejes
view(2);
title ('Placa');
xlabel 'EJE X'
ylabel 'EJE Y'
axis equal
hold off
3 Función Potencial del Fluido
Una vez definida la región que se va a estudiar, se analizará la velocidad de las partículas del fluido, que viene dada por su función potencial.
La función potencial establecida es [math] \varphi (\rho ,\theta)=(\rho +\frac{1}{\rho})\cos (\theta )+\sqrt{2}\log({\rho}) [/math]. Se representará la función potencial, para luego examinar la velocidad de las partículas y el movimiento del fluido.
rho=linspace(1,5,30); %Se define rho y theta
th=linspace(0,2*pi,30);
[U,V]=meshgrid(rho,th); %Se crea el mallado
hold on
X=U.*cos(V); %Se parametriza la superficie
Y=U.*sin(V);
f=@(rho,th)(rho+(1./rho)).*cos(th) +sqrt(2).*log(rho); %Se define la función potencial
Z=f(U,V); %Se aplica la función potencial
surf(X,Y,Z); %Se representa la función
plot(1*cos(th),1*sin(th),'k','lineWidth',3); %Se representa el obstáculo
view(2);
axis([-5,5,-5,5]);
colorbar;
title ('Función potencial');
xlabel ('EJE X');
ylabel ('EJE Y');
axis equal
hold off
4 Campo de Velocidades del Fluido
La velocidad de las partículas del fluido viene dada por el gradiente de la función potencial. Esto determina el valor y la dirección de máximo crecimiento de la función potencial.
El gradiente de una función en coordenadas cilíndricas se obtiene a partir de la siguiente fórmula: \(\nabla f(\rho, \theta, z) = \frac{\partial f}{\partial \rho} \vec e_\rho + \frac{1}{\rho} \frac{\partial f}{\partial \theta} \vec e_\theta + \frac{\partial f}{\partial z} \vec e_z\)
Siguiendo esta fórmula el gradiente de la función potencial determinada es [math] \vec u=\nabla \varphi=[(1-\frac{1}{\rho^2})\cos(\theta)+\frac{\sqrt{2}}{\rho}]\vec e_\rho-[(1+\frac{1}{\rho^2})\sin(\theta)]\vec e_\theta [/math]
Como la representación se realiza en coordenadas cartesianas, se debe pasar el gradiente de la función potencial a estas. Para ello usamos la matriz de cambio de base de cilíndricas a cartesianas y la multiplicamos por el vector dado en coordenadas cilíndricas.
Matriz de cambio de base:
\begin{bmatrix} \cos(\theta) & -\sin(\theta) & 0 \\ \sin(\theta) & \cos(\theta) & 0 \\ 0 & 0 & 1 \end{bmatrix}
Tras el cambio obtenemos;
[math] \vec u=\nabla \varphi=\left( (1 - \frac{1}{{\rho^2}}) \cdot \cos(\theta)^2 + \cos(\theta) \cdot \frac{\sqrt{2}}{\rho} - (1 + \frac{1}{{\rho^2}}) \cdot \sin(\theta)^2 \right) \vec i + \left( (1 - \frac{1}{{\rho^2}}) \cdot \sin(\theta) \cdot \cos(\theta) + \sin(\theta) \cdot \frac{\sqrt{2}}{\rho} + (1 + \frac{1}{{\rho^2}}) \cdot \sin(\theta) \cdot \cos(\theta) \right) \vec j [/math]
A continuación, se representa el campo de velocidades.
rho=linspace(1,5,30); %Definimos rho y theta
th=linspace(0,2*pi,30);
[U,V]=meshgrid(rho,th); %Creamos la malla
X=U.*cos(V); %Parametrizamos la superficie
Y=U.*sin(V);
f=@(rho,th)(rho+(1./rho)).*cos(th) + sqrt(2).*log(rho); %Definimos la función potencial
Z=f(U,V); %Aplicamos la función
contour(X,Y,Z,30); %Dibujamos las curvas de nivel
hold on
%Definimos las componentes X e Y del gradiente
Cx=(1-(1./U.^2)).*cos(V).^2 +cos(V).*(sqrt(2)./U)+ (1+(1./U.^2)).*sin(V).^2;
Cy=(1-(1./U.^2)).*sin(V).*cos(V)+sin(V).*(sqrt(2)./U)- (1+(1./U.^2)).*sin(V).*cos(V);
quiver(X,Y,Cx,Cy); %Dibujamos el campo de velocidades
plot(1*cos(th),1*sin(th),'k','lineWidth',1); %Representación del obstáculo
axis([-5,5,-5,5]);
colorbar; %Añadimos una barra de color
axis equal
view(2);
hold off
4.1 Interpretación y Valor de [math]\vec u [/math]
Utilizando la gráfica anterior, se puede visualizar como el campo de velocidades [math] \vec u=\nabla \varphi [/math], representado por las flechas, es ortogonal en todo punto a las curvas de nivel de [math] \varphi [/math].
Para comprobar matemáticamente esta observación, se elige un vector [math] \vec n [/math] cualquiera normal a los puntos del obstáculo, y se verifica que cumple la fórmula [math]\vec{u}\cdot \vec{n}=0 [/math]. Esto demuestra que, efectivamente los vectores son ortogonales entre sí.
Se elige un vector [math] \vec n [/math] ortogonal a las curvas de nivel, como puede ser [math]-\vec e_\theta [/math]. Se evalúa el gradiente [math] \vec u [/math] en el punto [math]\theta = 0[/math], y al multiplicarlo por el vector [math] \vec n [/math] elegido, se tiene que:
[math]-\vec e_\theta \cdot \left [(1-\frac{1}{\rho^2})\cos(\theta)+\frac{\sqrt{2}}{\rho} \right]\vec e_\rho = 0 [/math]
Se evalúa [math] \vec u [/math] en la frontera del obstáculo:
[math]\vec u(\rho,0,z)=\nabla \varphi(\rho,0,z)= [(1-\frac{1}{\rho^2})+\frac{\sqrt{2}}{\rho}]\vec e_\rho [/math]
4.2 Campo de velocidades lejos del obstáculo
4.3 Rotacional Nulo y Divergencia Nula
El rotacional y la divergencia del campo vectorial proporcionan información sobre las características físicas del fluido.
El operador rotacional indica la dirección y la velocidad del giro del campo en cada punto. Considera el movimiento del fluido y su tendencia a inducir rotación en cada punto.
[math]\nabla\times \vec{u}=\frac{1}{\rho }\begin{vmatrix} \vec e_\rho &\vec e_\theta &\vec e_z \\ \frac{\partial }{\partial \rho }& \frac{\partial }{\partial \theta }& \frac{\partial }{\partial z}\\ \cos \left (\theta \right ) \left ( 1-\frac{1}{\rho ^2} \right )& \rho \left [ -\frac{1}{\rho }sin \left ( \theta \right )\cdot\left ( \rho +\frac{1}{\rho } \right ) +\frac{\sqrt{2}}{\rho }\right ] & 0 \end{vmatrix}=[/math]
La divergencia en cambio, mide la tasa de expansión y contracción de un fluido en cada punto, calculando la diferencia entre el flujo entrante y saliente de la superficie.
5 Líneas de Corriente
Se va a proceder a calcular las líneas de corriente del campo u⃗, donde las líneas de corriente son aquellas que son tangentes al campo u⃗ en cada punto. Dichas líneas determinaran la trayectoria que sigue las partículas del fluido.
Para ello se calculará el campo que es ortogonal al campo u⃗. Este vector se asignará como v⃗ y se calculará con la siguiente fórmula: v⃗ = k⃗ x u⃗
En el apartado anterior se puede apreciar que el rotaciones es nulo y esto implica v⃗ sea irrotacional y además tiene un potencial escalar ψ cuyo gradiente es v⃗. Este potencial se conoce como función de corriente de u⃗.
Para calcular dicho potencial escalar ψ se resolverá este sistema de ecuaciones: ∇ψ= v⃗ , y se dibujarán las líneas ψ = cte.
