ANEXO I SOL BLASIUS

Vista previa del material en texto

AI - 1 
 
Anexo I.- Desarrollo y cálculo del modelo 
analítico. 
 
1.1.- Alcance del Anexo. 
En este anexo se desarrollara el modelo analítico de capa limite laminar utilizado como 
referencia para la comparación con nuestras soluciones numéricas. Además se expondrá la 
resolución de las mismas ecuaciones para nuestros casos particulares y las decisiones tomadas. 
Prácticamente la totalidad de dicha resolución se ha efectuado mediante el programa 
matemático MATLAB, en el que se ha implementado un código. Este código se reflejara 
prácticamente integro al final del anexo, y explicado por partes en cada apartado 
correspondiente. 
 
1.2.- Flujo laminar sobre placa plana. 
Para el problema de flujo laminar sobre placa plana empezaremos por considerar las siguientes 
hipótesis: 
 - Flujo estacionario 
 - Flujo incompresible 
 - Propiedades del fluido constantes 
 - Efecto de la viscosidad despreciable 
 - No existe gradiente de presiones a lo largo de la placa. 
Por tanto, tras asumir dichos puntos, las ecuaciones de la capa limite quedan reducidas a: 
Continuidad 
 
 
 
 
 
 
 
 
 
Cantidad de movimiento: 
 
 
 
 
 
 
 
 
 
 
 
 
Energía: 
 
 
 
 
 
 
 
 
 
 
 
 
Especies: 
 
 
 
 
 
 
 
 
 
 
 
 
Debemos notar a partir de estas 4 ecuaciones que el problema Hidrodinámico está desligado de 
la temperatura y de la concentración de especies, por tanto, competerá tan solamente a la 
ecuación 1 y 2. 
 
 
 
AI - 2 
 
1.2.1.- Problema hidrodinámico. 
 
Para la resolución de dicho problema utilizaremos el método de Blausius como se ha definido 
en la memoria. 
El primer paso será definir la función de flujo : 
 
 
 
 
 
 
 
 
Por tanto, con esto cumpliremos la ecuación de continuidad (1). Definimos dos nuevas 
variables, una dependiente y otra independiente: 
 
 
 √
 
 
 
 √
 
 
 
 
La definición de estas dos nuevas variables nos permitirá pasar en la ecuación de cantidad de 
movimiento, de una ecuación con derivadas parciales a una ecuación con derivadas ordinarias. 
 
La solución de Blausius es determinada como una solución de semejanza, y como una 
variable de semejanza. Se utiliza esta terminología porque, a pesar del crecimiento de la capa 
limite en función de nuestra posición en la placa, el perfil de velocidad permanece 
geométricamente semejante. Esta semejanza es de la forma: 
 
 
 (
 
 
) 
Donde es el espesor de la capa limite. Descubriremos de la solución de Blausius que el 
espesor varía como (
 
 
)
 
 
 , por tanto: 
 
 
 
Así, el perfil de velocidades vendrá únicamente determinado por la variable de semejanza la 
cual depende de x e y: 
 
 
 
 
 
 
 
 
 √
 
 
 
 
 √
 
 
 
 
 
 
 
 
 
y 
 
 
 
 ( √
 
 
 
 
 
 
 
 √
 
 
 ) 
 
 
 
 √
 
 
( 
 
 
 ) 
 
 
 
Ahora se deriva los componentes de la velocidad: 
 
 
 
 
 
 
 
 
 
 
 
 √
 
 
 
 
 
AI - 3 
 
 
 
 
 
 
 
 
 
 
 
Si estas tres últimas derivadas las sustituimos en la ecuación de conservación de la cantidad de 
movimiento obtenemos: 
 
 
 
 
 
 
 
De esta manera, el problema de capa limite hidrodinámica se ve reducido a un problema no 
lineal con derivadas ordinarias de tercer orden. Las condiciones de contorno para este problema 
serán: 
 
O bien: 
 
 
|
 
 
 
 
|
 
 
Estas ecuaciones serán resueltas por integración numérica como se mostrará en el capítulo 1.3. 
 
 
1.2.2.- Problema térmico. 
 
Así, una vez que resolvamos dicho problema, podremos resolver el térmico ya que este sí es 
dependiente de las velocidades halladas en el hidrodinámico. 
 
Una vez tenemos conocido el perfil de velocidades de la capa limite hidrodinámica 
resolveremos ahora la ecuación de la energía. 
 
Comenzamos con la introducción del concepto de temperatura adimensional 
 y asumiendo que la solución será del tipo 
 . 
Sustituyendo en la ecuación de la energía obtenemos: 
 
 
 
 
 
 
 
 
 
Debe notarse la dependencia que tiene la solución térmica con la solución hidrodinámica 
mediante la aparición de la variable f en la ecuación. Las condiciones de contorno pertinentes 
serán: 
 
AI - 4 
 
De nuevo, resolveremos la ecuación previa por integración numérica para un Prandtl 
determinado. 
 
 
 
|
 
 
 
 
 
 
1.3.- Resolución numérica. 
 
En este capítulo daremos resolución numérica a los problemas anteriormente planteados ya que 
estos serán necesarios para la comparación con los resultados obtenidos mediante el programa 
de CFD. 
 
Dicha resolución se ha implementado mediante la herramienta informática MATLAB. En 
primer lugar debemos mencionar que ambos problemas están compuestos por ecuaciones 
diferenciales, las cuales integraremos numéricamente mediante un método de Newton-Prince 
mejorado. Este ya se encuentra implementado como función en MATLAB como “ode45”. 
 
La principal contrariedad que se nos plantea con esta resolución es que dicha función solo 
trabaja con derivadas de primer orden, siendo las nuestras de segundo y tercero: 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
No obstante realizaremos el siguiente artificio para poder llevar a cabo la implementación de 
dicha ecuación. 
 
Esto consiste simplemente en transformar la ecuación de tercer orden a tres ecuaciones de 
primer orden de la manera mostrada en el código: 
 
 
function df = BlasiusFunc(eta,f) 
% f() es el vector de f y sus derivadas f = f(1) f' = f(2) f" = f(3) 
df = zeros(3,1); 
% df() son las derivadas de f() 
df(1) = f(2); 
df(2) = f(3); 
df(3) = -f(1)*f(3)/2; 
 
 
A su vez, esta será la función que implementaremos en “ode45”. 
Para el caso térmico se hará de manera análoga, pero al depender del hidrodinámico se 
conservara la primera parte de la función y por consiguiente la nomenclatura no será muy 
adecuada. Esto queda subsanado por las anotaciones en el propio código. 
 
 
function df = BlasiusFunc(eta,f) 
% f() es el vector de f y sus derivadas: f = f(1) f' = f(2) f" = f(3) 
df = zeros(3,1); 
AI - 5 
 
% df() son las derivadas de f() 
df(1) = f(2); 
df(2) = f(3); 
df(3) = -f(1)*f(3)/2; 
 
%Definimos el Prandlt para nuestro caso particular 
Pr=1006.43*0.002/1; 
%Subsanacion de la nomenclatura 
%df(4)=dT(1)=T* 
%df(5)=dT(2)=dT*/dy 
% dT(1)=T(2); 
% dT(2)=-(Pr/2)*f(1)*T(2); 
%f(4)=T(1) 
df(4)=f(5); 
df(5)=-(Pr/2)*f(1)*f(5); 
 
 
Debemos remarcar que ahora que nuestros problemas se han visto reducidos a 3 y 2 ecuaciones 
de primer orden respectivamente, debemos transformar también las condiciones de contorno de 
ambos problemas para adecuarlos a la propia función “ode45”. 
 
En primer lugar, para el problema hidrodinámico, debemos transformar la condición de 
contorno , en la condición inicial de y(3). Para ello generamos el siguiente código 
matlab en el que se representará todo el rango de puntos equivalentes entre estas dos variables y 
elegiremos el necesario (resolución inversa). 
 
 
for i=1:4000 q(i)=i/4000; 
[t,x]=ode45(@BlasiusFunc, [0,10],[0 0 q(i)]); %q(i) será la condición a establecer. 
r(i)=x(length(x),2); %La segunda derivada es la velocidad adimensional. 
 
end; 
figure(1); 
plot(q,r);grid;xlabel('Cond. Iniciales de y(3) equivalentes');ylabel('Valores de T(inf)'); 
 
AI - 6 
 
 
 
Una vez determinada esta transformación de la condición de contorno, estaremos en disposición 
de determinar también la necesaria para el caso térmico ( en la condición inicial de 
y(5) según se definió la función previamente). Por tanto, seguiremosunos pasos análogos a los 
anteriores para determinar esta: 
 
for i=1:4000 q(i)=i/4000; 
[t,x]=ode45(@BlasiusFunc, [0,10],[0 0 0.3322 0 q(i)]); %q(i) será la condición a establecer. 
r(i)=x(length(x),4); %Temp estabilizada 
end; 
figure(1); 
plot(q,r);grid;xlabel('Cond. Iniciales de y(3) equivalentes');ylabel('Valores de T(inf)'); 
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
0
0.5
1
1.5
2
2.5
X: 0.3322
Y: 1
Cond. Iniciales de y(3) equivalentes
V
a
lo
re
s
 d
e
 u
(i
n
f)
AI - 7 
 
 
Por tanto ahora ya tenemos determinada completamente ambos problemas de forma analítica de 
los cuales se muestran a continuación unas gráficas de los mismos. 
 
 
 
1.4.- Postratamiento de datos. 
El paso siguiente a la obtención de soluciones para los modelos analíticos y teniendo como 
datos los puntos experimentales a contrastar es el postratamiento de estos. Para ello se 
continuará con la herramienta MATLAB. 
 
1.4.1.- Implementación de datos CFD. 
La metodología seguida para la implementación de los datos obtenidos por CFD se describirá a 
continuación. 
 
En primer lugar, las soluciones se han exportado a formato Excel en forma de tablas 
[Velocidad/Temperatura; Posición eje coordenadas]. Para la implementación de estos datos en 
MATLAB seguiremos las siguientes tácticas en aras de crear un patrón general que facilite el 
0 0.2 0.4 0.6 0.8 1 1.2 1.4
0
5
10
15
20
25
df/dy (m/s)
e
ta
Perfil de velocidades analitico
0 0.2 0.4 0.6 0.8 1 1.2 1.4
0
5
10
15
20
25
T* (adim)
e
ta
Perfil de temperaturas analitico para Pr=2.3
AI - 8 
 
tratamiento de datos. Se definen funciones como las que siguen a continuación, en las que los 
valores se introducen como vectores y matrices y son llamados desde el programa principal. 
 
De esta manera evitamos sobrecargas en el programa principal ya que este no contiene 
directamente todos los valores, y además lo dota de mayor versatilidad al poder llamar a 
diferentes valores en cualquier momento. 
 
function [u,y]=disflu_A_I 
% 1,2,3 son el numero de muestras, u e y son sus componentes en la malla. 
u(1,:)=[]'; 
u(2,:)=[]'; 
u(3,:)=[]'; 
y(1,:)=[]'; 
y(2,:)=[]'; 
y(3,:)=[]'; 
 
 
Para el caso térmico seria de forma análoga: 
 
function [u,y,T,yT]=disflu_A_I 
% 1,2,3 son el numero de muestras, u e y son sus componentes en la malla. 
u(1,:)=[]'; 
u(2,:)=[]'; 
u(3,:)=[]'; 
y(1,:)=[]'; 
y(2,:)=[]'; 
y(3,:)=[]'; 
T(1,:)=[]'; 
T(2,:)=[]'; 
T(3,:)=[]'; 
yT(1,:)=[]'; 
yT(2,:)=[]'; 
yT(3,:)=[]'; 
 
 
 
1.4.2.- Contraste entre soluciones. 
A lo largo de todo el desarrollo de la metodología de obtención de la malla se requieren 
diferentes comparaciones entre soluciones (ver convergencia, errores relativos, etc.). 
 
Esto se hará principalmente de dos maneras. La primera se trata de jugar simplemente con la 
función “plot”, con lo que no entraremos en más detalle. 
Sin embargo, existe un segundo método orientado a las “Graficas de 45º” en las que el número 
de elementos a comparar en el eje X e Y deben tener en común el valor de Y en la malla. 
 
Debido a la heterogeneidad que tenemos entre el número de nodos entre las diferentes mallas 
esto no sería posible salvo por la solución propuesta. 
 
Procederemos así a realizar una interpolación con las valores “dato” para obtener un nuevo 
vector de velocidades o temperaturas en función de las Y’s convenientes. 
 
AI - 9 
 
Se ha creído conveniente emplear la función “interp1” en modo “Spline” como se muestra en el 
código a continuación. 
 
%siendo k, el número de la muestra, e y el nuevo vector de elementos de y que nos interesa 
 
x1=interp1(Yexp(k,:),Uexp(k,:),y,'spline'); 
x2=interp1(YexpII(k,:),UexpII(k,:),y,'spline');