Visualización de Campos Escalares y Vectoriales en Fluidos en Canales (Grupo D18)
Contenido
1 Introducción a la Teoría de la Dinámica de Fluidos
La teoría de la dinámica de fluidos ha sido un problema a resolver hasta las ecuaciones de Navier-Stokes. El primero en abordar este problema fue Newton, que aplico sus leyes de la mecánica a los movimientos de los fluidos, seguido de Euler, que consiguió definir una serie de ecuaciones diferenciales que rigen el movimiento de un fluido ideal. No fue hasta Navier y Stokes que consiguieron introducir la viscosidad y la temperatura en las ecuaciones definidas con anterioridad. Como resultado se obtuvo la siguiente ecuación:
[math]
\vec u \cdot \nabla \vec u + \nabla p = \mu \Delta \vec u,
[/math]
Siendo [math] \vec u (x,y) [/math] el campo vectorial de las velocidades de las partículas, p(x,y) el campo escalar de las presiones y [math] \mu [/math] es el coeficiente de rozamiento de las superficies de contacto entre el fluido y una superficie.
Para que el lector comprenda mejor las ecuaciones las aplicaremos a un ejemplo práctico en R² explicando sus consecuencias. Este se especifica a continuación:
CAMPO DE VELOCIDADES: [math] \vec u(x,y)= y(1-y)(p_1-p_2)/(2\mu) \vec i [/math]
CAMPO DE PRESIONES: [math] p(x,y)=p_1 + (p_2-p_1)(x-1) [/math]
En un canal plano de dimensiones [0,4]x[0,1]
2 Razonamientos con los Campos Escalares y Vectoriales
Partimos del cálculo del gradiente del campo de velocidades y lo multiplicamos por el campo, al ser un producto escalar el orden de los factores no altera el resultado final.
[math]
\vec u \cdot \nabla \vec u + \nabla p = \mu Δ \vec u
[/math]
Calculamos ambos miembros de la ecuación:
[math]
\nabla \vec u = \frac{ \partial u_{i}}{ \partial x^{j} } = \begin{bmatrix} \frac{\partial u_1}{\partial x} & \frac{\partial u_1}{\partial y} \\\frac{\partial u_2}{\partial x} & \frac{\partial u_2}{\partial y} \end{bmatrix} =\begin{bmatrix} 0 & \frac{(1-2y)(p_{1}-p_{2})}{(2\mu)} \\ 0 & 0 \end{bmatrix}
[/math]
[math]
\vec u \cdot \nabla \vec u = \nabla \vec u \cdot \vec u = \begin{bmatrix} \frac{\partial u_1}{\partial x} & \frac{\partial u_1}{\partial y} \\\frac{\partial u_2}{\partial x} & \frac{\partial u_2}{\partial y} \end{bmatrix} \begin{bmatrix} u_1 \\ u_2 \end{bmatrix} = \begin{bmatrix} 0 & \frac{(1-2y)(p_{1}-p_{2})}{(2\mu)} \\ 0 & 0 \end{bmatrix} \begin{bmatrix} y(1-y)(p_1-p_2)/(2\mu) \\ 0 \end{bmatrix} = \vec 0
[/math]
En el primer miembro solo queda el gradiente del campo escalar de presiones:
[math] \nabla p = \frac{ \partial p }{ \partial x } \vec i + \frac{ \partial p }{ \partial y} \vec j = (p_{2}-p_{1}) \vec i [/math] Al otro lado de la igualdad calculamos el laplaciano de [math]\vec u[/math]: [math] Δ \vec u = \nabla \nabla \vec u = \frac{ \partial^2 \vec u }{ \partial x^2 } + \frac{ \partial^2 \vec u }{ \partial y^2 }= \frac{ \partial^2 u_1 }{ \partial x^2 } \vec i + \frac{ \partial^2 u_2 }{ \partial x^2 } \vec j + \frac{ \partial^2 u_1 }{ \partial y^2 } \vec i + \frac{ \partial^2 u_2 }{ \partial y^2 } \vec j = \frac { p_2-p_1}{\mu} \vec i [/math]
Entonces al multiplicar por el coeficiente de rozamiento simplificamos y resulta
Quedando demostrada entonces las Ecuaciones de Navier-Stokes
2.2 Condión de incompresibilidad
Por otro lado la condición de incompresibilidad de un fluido viene dada por la igualdad [math] \nabla \cdot \vec u =0 [/math], donde se tiene que [math] \nabla \cdot \vec u = Tr(\nabla \vec u) [/math]
[math]
Tr(\nabla \vec u)= Tr \begin{bmatrix} \frac{\partial u_1}{\partial x} & \frac{\partial u_1}{\partial y} \\\frac{\partial u_2}{\partial x} & \frac{\partial u_2}{\partial y} \end{bmatrix} = Tr \begin{bmatrix} 0 & \frac{(1-2y)(p_{1}-p_{2})}{(2\mu)} \\ 0 & 0 \end{bmatrix} = 0
[/math]
3 Representación de los Campos escalares y Vectoriales
3.1 Campos de Velocidad y Presión
La visualización de los campos resulta muy útil para estudiar el comportamiento de los fluidos. Por ello se dibujarán el campo de velocidad [math] \vec u(x,y)[/math] y el campo de presiones [math] \vec p(x,y)[/math] , tomando los valores [math]p_1=2[/math], [math]p_2=1[/math] y [math]\mu=1[/math]
CAMPO DE VELOCIDADES:[math] \vec u(x,y)= y(1-y)/2 \vec i [/math] CAMPO DE PRESIONES:[math] p(x,y)= 3-x[/math]
Se observa en las visualizaciones, que en los extremos del canal la velocidad es nula, debido al rozamiento con las paredes del canal que se ha tenido en cuenta, mientras que en el medio la velocidad es máxima.
También se ve que la presion disminuye a media que X aumenta.
3.2 Lineas de Corriente
Las lineas de Corriente representa la dirección en la que se desplaza el fluido, siendo esta de gran interés. Estas lineas se obtienen de forma matemática como el campo escalar que es potencial del campo ortogonal al campo de velocidades, obteniendose las lineas tangentes a los vectores de [math] \vec u(x,y)[/math] en cada punto.
Llamemos [math] \vec v(x,y)[/math] al campo ortogonal a [math] \vec u(x,y)[/math].
[math] \vec v(x,y)= y(1-y)/2 \vec j [/math]
Luego el campo Escalar Potencial de [math] \vec v(x,y)[/math] es [math] \psi [/math], donde se cumple que
[math]
\nabla \psi = \vec v(x,y)= y(1-y)/2 \vec j
[/math]
[math]
\frac{\partial \psi}{\partial x}=0 \vec i
[/math]
[math]\frac{\partial \psi}{\partial y}= y(1-y)/2 \vec j
[/math]
[math]\psi = \int (y(1-y)/2) \partial y
[/math]
[math]\psi = \frac {(3y^2 -2y^3)}{2}
[/math]
Al visualizarlo queda como muestra la imagen, que son un conjunto de lineas paralelas entre sí y al eje OX.
3.3 Rotacional
Un gráfico que muestre el Rotacional del campo de Velocidades, nos muestra la velocidad a la que las partículas giran sobre si mismas debido a las diferencias de velocidad entre las que están a su alrededor. En el gráfico de la imagen se muestra el módulo del rotacional de [math] \vec u(x,y)[/math].
[math]
\nabla\times \vec u = \left|
\begin{matrix}
\hat x & \hat y & \hat z \\
& & \\
\frac{\partial}{\partial x} & \frac{\partial}{\partial y} & \frac{\partial}{\partial z}
\\ & & \\
F_x & F_y & F_z
\end{matrix}\right|= \left|
\begin{matrix}
\hat x & \hat y & \hat z \\
& & \\
\frac{\partial}{\partial x} & \frac{\partial}{\partial y} & \frac{\partial}{\partial z}
\\ & & \\
\frac{y(1-y)}{2} & 0 & 0
\end{matrix}\right|
[/math]
[math]
\nabla\times \vec u = (\frac{1}{2}-y) \vec k
[/math]
[math]
\left|
\nabla\times \vec u
\right| = \left|
(\frac{1}{2}-y) \vec k
\right| = \left|
(\frac{1}{2}-y)
\right|
[/math]
Los puntos de mayor rotacional, donde el giro de la particula es mayor, es en los bordes del canal, causado sobre todo por el rozamiento del canal y el avance de las particulas de alrededor. Generado por una gran diferencia que hace que el rotacional sea mayor que en el resto del canal, se repite a lo largo de él ya que al igual que en otros campos y gráficas, presenta una simetria axial respecto al eje del canal.
3.4 Campo de Temperaturas
Si la temperatura del fluido del canal viene dada por la función [math]T(x,y)=(x-1)^2-y^2[/math] al reprentarlo se obtiene el gráfico indicado. También se ha representado el gradiente, donde se observa que los vectores son ortogonales a las lineas isotermas.
[math]
\nabla T = \frac{ \partial T }{ \partial x } \vec i + \frac{ \partial T }{ \partial y } \vec j = 2(x-1) \vec i - 2y \vec j
[/math]
4 Anexo:programas utilizados en Matlab
El programa Matlab ha sido utilizado para la representación de todas las gráficas expuestas en este artículo. A continuación se muestran todos lo programas.
4.1 Mallado del Canal en un rango [0,4]x[-1,2]
h=0.1;
u=0:h:4;
v=0:h:1; %Dividimos en porciones de medida h=0.1 las dimensiones del canal
[uu,vv]=meshgrid(u,v); %Creamos las matrices de la red
figure(1) %Creamos la primera imagen
mesh(uu,vv,0*vv) %creamos una red en la zona del canal, indicando una
%medida para la coordenada Z en función alguna variable
axis([0,4,-1,2]) %Delimitamos el gráfico mostrado
view(2) %Pedimos que nos muestre el gráfico
print('malla.png') %Pedimos que cree una imagen del gráfico mostrado en formato .png
4.2 Representación del Campo Vectorial de Velocidades
x=linspace(0,4,8); %Dividimos en 8 porciones la dimensión x del canal
y=linspace(0,1,11); %Dividimos en 11 porciones la dimensión Y del canal
%El número de divisiones divisiones que se han puesto, son los más
%indicados para una mejor visualización
[X,Y]=meshgrid(x,y); %Creamos las matrices de la red
U=(Y.*(1-Y))/2; %determinmos la función de los vectores i
V=0*U; %y la de los vectores j, en función de una variable
quiver(X,Y,U,V) %Pedimos que nos represente el campo de vectores
axis([0,4,-1,2]) %Delimitamos el gráfico mostrado
print('velocidades.png') %Pedimos que cree una imagen del gráfico mostrado en formato .png
4.3 Representación del Campo Escalar de Presiones
h=0.1;
u=0:h:4;
v=0:h:1; %Dividimos en porciones de medida h=0.1 las dimensiones del canal
[uu,vv]=meshgrid(u,v); %Creamos las matrices de la red
Z=3-uu; %Indicamos la función de la presión
h=pcolor(uu,vv,Z); %Pedimos que nos represente el campo escalar mediante una gama de
%colores
axis([0,4,-1,2]) %Delimitamos el gráfico mostrado
print('presiones.png') %Pedimos que cree una imagen del gráfico mostrado en formato .png
4.4 Lineas de corriente
h=0.1;
u=0:h:4;
v=0:h:1; %Dividimos en porciones de medida h=0.1 las dimensiones del canal
[uu,vv]=meshgrid(u,v); %Creamos las matrices de la red
Z=(((vv.^2)./4)-((vv.^3)./6)); %Indicamos la función de psi de las lineas de corriente
contour(uu,vv,Z) %Pedimos que nos represente las lineas de nivel del campo
%escalar
axis([0,4,-1,2]) %Delimitamos el gráfico mostrado
print('corriente.png') %Pedimos que cree una imagen del gráfico mostrado en formato
%.png
4.5 Gráfico del módulo del Rotacional
x=linspace(0,4,8); %Dividimos en 8 porciones la dimensión x del canal
y=linspace(0,1,21); %Dividimos en 21 porciones la dimensión 4 del canal
[uu,vv]=meshgrid(x,y); %Creamos las matrices de la red
Z=abs(0.5-vv); %Indicamos la función del módulo del rotacional
h=pcolor(uu,vv,Z);; %Pedimos que nos represente el campo escalar mediante una gama de
%colores
axis([0,4,-1,2]) %Delimitamos el gráfico mostrado
print('rotacional.png')%Pedimos que cree una imagen del gráfico mostrado en formato .png
4.6 Lineas Isotermas y Campo Gradiente de Termperaturas
h=0.1;
u=0:h:4;
v=0:h:1; %Dividimos en porciones de medida h=0.1 las dimensiones del canal
[X,Y]=meshgrid(x,y); %Creamos las matrices de la red
T=((X-1).^2)-Y.^2; %Indicamos la función de la temperatura
U=2.*(X-1); %la función del gradiente en el vector i
V=-2.*Y; %la función del gradiente en el vector j
contour(X,Y,T) %Pedimos que nos represente las lineas de nivel del campo
%escalar de temperaturas
hold on
quiver(X,Y,U,V) %junto a las lineas de nivel pedimos que represente el campo
%vectorial gradiente del campo de temperaturas.
axis([0,4,-1,2]) %Delimitamos el gráfico mostrado
print('rotacional.png') %Pedimos que cree una imagen del gráfico mostrado en formato
%.png
hold off
