Flujo de Couette entre dos tubos concéntricos (Grupo 25)

Trabajo realizado por estudiantes
Título	Flujo de Couette entre dos tubos concéntricos . Grupo 25
Asignatura	Teoría de Campos
Curso	2025-26
Autores	Andrés Gaitan, Pablo Casado, Esteban Gabaldon, Carlos Ascanio, Aarón Perez
Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura

Flujo de Couette entre dos tubos concéntricos

Introducción

El flujo de Couette entre dos tubos concentricos, designa el movimiento de un fluido viscoso entre dos cilindros coaxiales de diferente radio, siendo dominado así el movimiento por la rotación de los cilindros alrederor de una eje común. En ingenieria, el flujo de Couette es sumamente importante ya que es uno de los pocos fluidos viscosos de los que se conoce la solución exacta de la ecuacion de Navier-Stokes, volviendolo fundamental en el estudio de propiedades de los fluidos y estabilidad de los flujos.

En este trabajo partiremos de un flujo generado por la rotación en sentido antihorario del cilindro interior con velocidad angular constante [math]\vec{w}_i[/math] y la rotación en sentido horario del cilindro exterior con velocidad angular [math]\vec{w}_e[/math]. Ambos cilindros tienen su eje [math]OX_3[/math] y sus radios [math]\rho[/math]=2 para el cilindro exterior y [math]\rho[/math]=1 para el interior. A partir de ello, con el campo de velocidades entre ambos cilindros, estudiaremos el flujo de temperaturas resultante, el gradiente y el caudal, entre otras.

1 Mallado de la sección transversal

En primer lugar, se toma la sección de ambos cilindros en el plano [math]x_3[/math]=0 y se representa en matlab de acuerdo al parámetro angular [math]\theta \in [0,2\pi][/math] y a los parámetros cartesianos [math](x,y) \in [-2,2] \times [-2,2][/math].

% MALLADO DE LA SECCIÓN TRANSVERSAL
figure; 
rho=1:0.1:2;                                  % SE ESTABLECEN PARAMETROS ALGULARES Y CARTESIANOS                    
theta_malla = linspace(0, 2*pi, 50);             
[RHO,THETA]=meshgrid(rho,theta_malla);        % SE OBTIENE LA MALLA         
rho_int=1;
rho_ext=2;
hold on
grid on
x=RHO.*cos(THETA);                               
y=RHO.*sin(THETA);                               
mesh(x,y,0*y,'EdgeColor','g');                % SE REPRESENTA LA MALLA
theta_suave = theta_malla; 
x_ext = rho_ext.*cos(theta_suave);            % SE OBTIENEN LOS CILINDROS INTERIOR Y EXTERIOR EN PLANTA
y_ext = rho_ext.*sin(theta_suave);
x_int = rho_int.*cos(theta_suave); 
y_int = rho_int.*sin(theta_suave);
axis([-2, 2, -2, 2]);                         % SE LIMITAL LOS EJES                  
axis equal;                                      
plot(x_ext, y_ext, 'b', 'LineWidth', 3);      % SE REPRESENTAN AMBOS CILINDROS  
plot(x_int, y_int, 'b', 'LineWidth', 3);         
view(2)                                       % MALLADO DE LA SECCIÓN TRANSVERSAL - MALLA CIRCULAR Y BAJA DENSIDAD
hold off

2 Calculo de velocidades. Ecuación de Navier-Stokes

La velocidad de las partículas del fluido viene dada por la expresión [math]\vec{u}(\rho, \theta) = f(\rho)\vec{e}_{\theta}[/math], donde la presión p es constante. Sabiendo que (\vec{u},p) cumple la ecuacion de Navier-Stokes estacionaria se obtiene:

[math](\vec{u} \cdot \nabla)\vec{u} + \nabla p = \mu \Delta\vec{u}[/math]

Siendo \mu el coeficiente de viscosidad del fluido.

Además, como la presión p es constante, su gradiente sera nulo: [math]\nabla p = 0[/math], y si aparte despreciamos el primer término de la fórmula anterior obtenemos que:

[math]µ∆\vec{u} =\vec{0}[/math]

Así, se llega al Laplaciano de un campo vectorial que se desarrollara en la base cilíndrica de acuerdo a la siguiente expresión:

[math]\Delta\vec{u} = \nabla(\nabla \cdot \vec{u}) - \nabla \times (\nabla \times \vec{u})[/math]

Se desarrolla el primer termino, que es el gradiente de una divergencia en coordenadas cilindricas:

[math]\nabla \cdot \vec{u} = \frac{1}{\rho} \frac{\partial}{\partial \rho}(\rho u_{\rho}) + \frac{1}{\rho} \frac{\partial u_{\theta}}{\partial \theta} + \frac{\partial u_{z}}{\partial z}[/math]

Si se particulariza al ejercicio se obtiene que:

[math]\nabla \cdot \vec{u} = \frac{1}{\rho} \left[ \frac{\partial}{\partial \rho}(0) + \frac{\partial}{\partial \theta}(f(\rho)) + \frac{\partial}{\partial z}(0) \right] = 0[/math]

Llegandose asi a la conclusión de que el primer término es igual al vector nulo: [math]\nabla(\nabla \cdot \vec{u}) = \vec{0}[/math]

Ahora se desarrolla el segundo termino, el rotacional del rotacional de [math]\vec{u}[/math], para ello se empieza en primer lugar con el primer rotacional:

[math]\nabla \times \vec{F}(\rho, \theta, z) = \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_{\rho} & \rho F_{\theta} & F_{z} \end{vmatrix}[/math]

Así, particularizando al ejercicio se obtiene:

[math]\nabla \times \vec{u} = \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} = \left( \frac{1}{\rho} \frac{\partial}{\partial \rho}(\rho f(\rho)) \right) \vec{e}_{z}[/math]

Tras ello, se desarrolla el siguiente rotacional:

[math]\nabla \times (\nabla \times \vec{u}) = - \frac{1}{\rho} \det\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 & \frac{1}{\rho} \frac{\partial}{\partial \rho} (\rho f(\rho)) & 0\end{vmatrix} = -\left[ \frac{\partial}{\partial \rho}\left(\frac{1}{\rho}\frac{\partial}{\partial \rho}(\rho f(\rho))\right) \right] \vec{e}_{\theta} = -\left[ \frac{\rho^2 \frac{\partial^2 f}{\partial \rho^2} + \rho \frac{\partial f}{\partial \rho} - f(\rho)}{\rho^2} \right] \vec{e}_{\theta}= -\left[ \frac{\partial^2 f}{\partial \rho^2} + \frac{1}{\rho} \frac{\partial f}{\partial \rho} - \frac{f(\rho)}{\rho^2} \right] \vec{e}_{\theta}[/math]

Por tanto, tras desarrollar todos los terminos se llega a que:

[math]\Delta\vec{u} = \left[ \frac{\partial^2 f}{\partial \rho^2} + \frac{1}{\rho} \frac{\partial f}{\partial \rho} - \frac{f(\rho)}{\rho^2} \right] \vec{e}_{\theta}[/math]

Ya, una vez conocidos todos los términos de la ecuación de Navier-Stokes, se puede resolver:

[math](\vec{u} \cdot \nabla)\vec{u} + \nabla p = \mu \Delta \vec{u} \Rightarrow 0 + \vec{0} - \mu \Delta \vec{u} = \vec{0} \Rightarrow \Delta\vec{u} = \vec{0}[/math]

De lo que se deduce la siguiente ecuación diferencial:

[math]\left[ \frac{\partial^2 f}{\partial \rho^2} + \frac{1}{\rho} \frac{\partial f}{\partial \rho} - \frac{f(\rho)}{\rho^2} \right] \vec{e}_{\theta} = 0[/math]

3 Representación del campo de velocidades

Partiendo de las siguientes expresiones; [math]\vec{u}(\rho,\theta) = f(\rho)\,\vec{e}_\theta,\quad f(\rho) = a\rho + \frac{b}{\rho}[/math]; y los valores anteriormente calculados; [math]a = \frac{2}{3}\omega \quad b = -\frac{2}{3}\omega[/math], si suponemos que [math]|\vec{\omega}_e| = |\vec{\omega}_i| = 1 \text{ y } \mu = 1[/math], llegamos a nuestro campo de velocidades:

[math]\vec u(\rho, \theta)=f(\rho)\vec e_\theta=\left ( \frac{4}{3}\rho -\frac{4}{3\rho}\right )\vec e_\theta[/math]