Flujo de Poiseuille Grupo 18-B

Trabajo realizado por estudiantes
Título	Flujo de Poiseuille (Grupo 18-B)
Asignatura	Alejandro Ramos García Alberto Monge Alejandro López Louciana Contreras
Curso	{{{3}}}
Autores	{{{4}}}
Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura

Vamos a considerar el flujo de un fluido incompresible en una tubería cilíndrica de radio 2 que supondremos centrada en el eje OZ, de la cual conocemos que la velocidad [math]\vec{u}(\rho,\theta,z)[/math] de las partículas viene dada por:

[math]\vec{u}(\rho,\theta,z)= f\left(\rho\right)\vec{e_{z}},[/math]

y su presión [math]p\left(x,y\right)[/math] que viene dada por:

[math]p\left(x,y\right)=p_{1}+\left (p_{2}-p_{1}\right)(z-1).[/math]

Posteriormente, representaremos el campo de velocidades y presiones gráficamente. Además, analizaremos los vectores ortogonales a [math]\vec{u}[/math], los puntos donde este es máximo, su rotacional y calcularemos el caudal que pasa por una de sus secciones transversales.

También, analizaremos una segunda cantidad física; la temperatura, la cual viene dada por el campo:

[math]T\left (\rho,\theta,z\right)=1 +\left (\rho-\frac{1}{2}\right)^{2}e^{-\left(z-1\right)^{2}},[/math]

de la cual, dibujaremos el campo de temperaturas y las curvas de nivel del campo. Luego, calcularemos su gradiente y analizaremos sus vectores gráficamente.

1 Mallado. Sección transversal de la tubería.

Dibujaremos un mallado, en dimensión 2, de la sección longitudinal de la tubería [math] x_{1} = 0 [/math]. Fijaremos los ejes en la región [math] \left ( \rho,z \right )\epsilon \left [ 0,3 \right ]\times \left [ 0,10 \right ]. [/math]

Sección transversal

clc;
clear all;
x=0:0.3:3;
y=0:0.3:10;
[xx,yy]=meshgrid(x,y)
figure(1)
mesh(xx,yy,0*xx)
axis([-2,5,0,10])
view(2);

2 Velocidad y presión

2.1 Ecuación de Navier-Stokes. Divergencia nula.

2.1.1 Ecuación de Navier-Stokes.

[math] \frac{1}{\rho} \frac{\partial }{\partial \rho }\left ( \rho \frac{\partial f\left ( \rho \right )}{\partial \rho } \right )=\frac{p_{2}-p_{1}}{\mu}\rightarrow \frac{\partial }{\partial \rho}\left ( \rho \frac{\partial f\left ( \rho \right )}{\partial \rho } \right )= \rho\frac{\left ( p_{2}-p_{1} \right )}{\mu} [/math]

[math] \int \frac{\partial }{\partial \rho}\left ( \rho \frac{\partial f\left ( \rho \right )}{\partial \rho } \right ) d\rho = \frac{\left ( p_{2}-p_{1} \right )}{\mu}\int\rho d\rho [/math]

[math] \rho \frac{\partial f\left ( \rho \right )}{\partial \rho } d\rho = \frac{\left ( p_{2}-p_{1} \right )}{\mu} \cdot \frac{\rho^{2}}{2} [/math]

[math] \rho\cdot{f}'\left(\rho\right)=\frac{p_{2}-p_{1}}{\mu}\cdot \frac{\rho^{2}}{2} [/math]

[math]{f}'\left(\rho\right)=\frac{p_{2}-p_{1}}{\mu}\cdot \frac{\rho}{2}[/math]

[math] \int \left ( \rho \frac{\partial f\left ( \rho \right )}{\partial \rho } \right ) d\rho = \frac{\left ( p_{2}-p_{1} \right )}{2\mu}\int\rho d\rho [/math] [math] f\left(\rho \right )=\frac{\left ( p_{2}-p_{1} \right )}{4\mu}\rho^{2} [/math]

2.1.2 Divergencia nula.

La divergencia de una campo vectorial en un punto, más concretamente, de un campo de velocidades de partículas, se puede interpretar como el cambio en la densidad en ese punto.

Al ser el agua un fluido incompresible, la divergencia debería ser nula en cualquier punto. Matemáticamente:

[math]\vec{u}(\rho,\theta,z)= f\left(\rho\right)\vec{e_{z}}[/math]

[math]\triangledown\cdot\vec{u} \left(\rho,\theta,z\right)=\frac{1}{\rho} \left ( \frac{\partial }{\partial \rho}\left ( \rho u_{\rho } \right )+\frac{\partial }{\partial \theta}\left ( u_{\theta } \right ) + \frac{\partial }{\partial z}\left (\rho u_{z}\right)\right) =\frac{1}{\rho}\left ( \frac{\partial }{\partial z}\left(\rho f\left(\rho\right)\right)\right)=0 [/math] (para cualquier punto).

2.2 Campo de presiones y campo de velocidades.

Considerando: [math] p_{1}=2, p_{2}=1 [/math] y [math]\mu=1, [/math] calcularemos el campo de presiones y el campo de velocidades.

El campo de presiones nos queda:

[math]p\left ( x,y \right )=p_{1}+\left ( p_{2}-p_{1} \right )\left ( z-1 \right )=2+\left ( 1-2 \right )\left ( z-1 \right )=3-z.[/math]

El campo de velocidades nos queda:

[math]\vec{u}(\rho,\theta,z)= f\left(\rho\right)\vec{e_{z}}=\frac{p_{2}-p_{1}}{4\mu}\rho ^{2}\vec{e_{z}}=-\frac{1}{4}\rho ^{2}\vec{e_{z}}[/math]

2.2.1 Representación del campo de presiones.

clc
clear all
%La función de presiones gráficamente es una recta.
z=0:0.1:10;
f=3-z;
plot(z,f)
xlabel('Variación de altura');
ylabel('Variación de presión');
title('Campo de presiones del fluido');

2.2.2 Representación del campo de velocidades.

clear;
clc;
x=0:0.1:3;
y=0:0.1:10;
[xx,yy]=meshgrid(x,y);
ux=-(1/4)*xx.^2;           %debido al cambio a coordenadas cartesianas
uy=-(1/4)*yy.^2;           %ídem
hold on
quiver(xx,yy,ux,uy)
axis([0,5,0,12])
hold off
view(2)
title('CAMPO DE VELOCIDADES DEL FLUIDO')

2.3 Líneas de corriente del campo.

Calcularemos el campo [math] \vec{v} [/math], el cual es ortogonal a [math]\vec{u}[/math] en cada punto.

[math] \vec{v}=\vec{e_{\theta}}\times \vec{u}=\begin{vmatrix} \vec{e_{\rho}} & \vec{e_{\theta}} & \vec{e_{z}} \\ 0 & 1 & 0\\ 1 & 0 & f\left ( \rho \right ) \end{vmatrix} = f\left ( \rho \right )\vec{e_{\rho }}=\frac{p_{2}-p_{1}}{4\mu}\rho ^{2}\vec{e_{\rho}} [/math]

El vector [math] \vec{v} [/math] es irrotacional al tener [math] \vec{u} [/math] divergencia nula. Matemáticamente se puede comprobar:

[math] \triangledown \times \vec{v}=\frac{1}{\rho}\begin{vmatrix} \vec{e_{\rho}}& \rho\vec{e_{\theta}} & \vec{e_{z}}\\ \frac{\partial }{\partial \rho}& \frac{\partial }{\partial \theta} & \frac{\partial }{\partial z}\\ u_{\rho} & \rho u_{\theta} & u_{z} \end{vmatrix} = \frac{1}{\rho}\begin{vmatrix} \vec{e_{\rho}}& \rho\vec{e_{\theta}} & \vec{e_{z}}\\ \frac{\partial }{\partial \rho}& \frac{\partial }{\partial \theta} & \frac{\partial }{\partial z}\\ f\left(\rho\right) & 0 & 0 \end{vmatrix} =\frac{1}{\rho}\left ( \frac{\partial f\left ( \rho \right ) }{\partial z} \rho \vec{e_{\theta}}- \frac{\partial f\left ( \rho \right ) }{\partial \theta}\vec{e_{z}}\right )=0 [/math]

Ahora calculamos [math]\psi[/math], función corriente de [math] \vec{u}. \left (\bigtriangledown\psi =\vec{v}\right) [/math]

[math] \triangledown \psi=\frac{\partial \psi}{\partial \rho} \vec{e_{\rho}} + \frac{\partial \psi}{\partial \theta} \frac{1}{\rho}\vec{e_{\theta}}+\frac{\partial \psi}{\partial z} \vec{e_{z}}=\frac{p_{2}-p_{1}}{4\mu}\rho ^{2}\vec{e_{\rho}} [/math]

[math] \psi=\frac{p_{2}-p_{1}}{4\mu}\int\rho ^{2}d\rho=\frac{p_{2}-p_{1}}{12\mu}\rho ^{3}=-\frac{1}{12}\rho^{3} [/math]

2.4 Puntos con velocidad máxima.

%Módulo del campo de velocidades
clc
clear all
rho=0:0.1:3;
f=(1/4)*(rho.^2);
plot(rho,f);
%El módulo es máximo cuando el rho es máximo, es decir, en rho=3
xlabel('aumento del rho')
ylabel('aumento de la velocidad')
title('Módulo del campo de velocidades');

2.5 Rotacional.

Podemos considerar que la interpretación del rotacional es el giro alrededor de un punto del fluido. Si el fluido fuera perfecto, el rotacional sería nulo, ya que, todas las velocidades serían iguales. En este caso, estamos trabajando con un fluido no perfecto, ya que, tiene coeficiente de viscosidad [math] \mu [/math], entonces el rotacional será distinto de cero.

Procedemos a calcularlo:

[math]\triangledown \times \vec{v}=\frac{1}{\rho}\begin{vmatrix} \vec{e_{\rho}}& \rho\vec{e_{\theta}} & \vec{e_{z}}\\ \frac{\partial }{\partial \rho}& \frac{\partial }{\partial \theta} & \frac{\partial }{\partial z}\\ u_{\rho} & \rho u_{\theta} & u_{z} \end{vmatrix} = \frac{1}{\rho}\begin{vmatrix} \vec{e_{\rho}}& \rho\vec{e_{\theta}} & \vec{e_{z}}\\ \frac{\partial }{\partial \rho}& \frac{\partial }{\partial \theta} & \frac{\partial }{\partial z}\\ 0 & 0 & f\left(\rho\right)\end{vmatrix} = -{f}'\left ( \rho \right )\vec{e_{\theta}} [/math]

Sustituímos con los valores de [math] p_{1}, [/math] [math] p_{2} [/math] y [math] \mu.[/math]

[math] \frac{p_{2}-p_{1}}{2\mu} = -\frac{\rho}{2} [/math]

Con esto podemos concluir que los puntos con mayor rotacional serán los puntos en los que la función tenga mayor pendiente.

2.5.1 Representación del campo del rotacional

clear;
clc;
x=0:0.1:3;
y=0:0.1:10;
[xx,yy]=meshgrid(x,y);
rot=abs(-(1/4)*xx.^2-(1/4)*yy.^2);
surf(xx,yy,rot)
colorbar
view(2)
axis([0,5,0,12])

2.6 Caudal por una sección transversal

Podemos calcular el caudal de una sección transversal integrando el campo mediante una integral (doble) de superficie.

[math] \int_{S}^{}\vec{u}\cdot d\vec{S} [/math]

3 Temperatura

3.1 Campo de temperaturas y curvas de nivel.

3.1.1 Representación del campo de temperaturas.

y=0:0.05:10;
z=-2:0.05:11;
[Y,Z]=meshgrid(y,z);
figure (1) 
p=1+(((Y.^2+Z.^2-(1/2)).^2).*exp(-((Z-1).^2)));
pcolor(Y,Z,p);
shading flat
grid on
axis([0,8,-1,10]);
colorbar

3.1.2 Representación de las curvas de nivel del campo de temperaturas.

y=0:0.05:8;
z=-2:0.05:10;
[Y,Z]=meshgrid(y,z);
p=1+(((Y.^2+Z.^2-(1/2)).^2).*exp(-((Z-1).^2)));
pcolor(Y,Z,p);
shading flat
grid on
axis([0,8,-1,10]);
colorbar 
contour(Y,Z,p,10,'k');

3.2 Gradiente.

El gradiente de la función de temperatura será:

[math] \bigtriangledown T\left ( \rho, \theta, z \right )=\frac{1}{4}e^{-\left ( z-1 \right )^{2}}\vec{e_\rho}+\left ( \rho-\frac{1}{2} \right )^{2}\cdot e^{-\left ( z-1 \right )^{2}}\cdot -2\left ( z-1 \right )\vec{e_z} [/math]

3.2.1 Representación del gradiente.

y=0:0.05:8;
z=0:0.05:1;
[Y,Z]=meshgrid(y,z);
 figure (1)
 p=1+((Z.^2).*exp(-((Y.^2)+(Z.^2)-1/2).^2));
[pY,pZ]=gradient(p);
hold on
quiver(Y,Z,pY,pZ)
axis([0,8,-1,2]);
shading flat
grid on
hold off

3.2.2 Representación curvas de nivel del gradiente.

y=0:0.05:8;
z=0:0.05:10;
[Y,Z]=meshgrid(y,z);
figure (1)
p=1+(((Y.^2+Z.^2-(1/2)).^2).*exp(-((Z-1).^2)));
[pY,pZ]=gradient(p);
hold on
quiver(Y,Z,pY,pZ)
contour(Y,Z,p,'k');
axis([0,8,-1,10]);
shading flat
grid on
hold off