Producción masiva de microburbujas monodispersas para aplicaciones reales evitando la microfluídica

Proyecto Fin de Carrera Ingeniería de Telecomunicación Formato de Publicación de la Escuela Técnica Superior de Ingeniería Autor: F. Javier Payán Somet Tutor: Juan José Murillo Fuentes Dep. Teoría de la Señal y Comunicaciones Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2013 Trabajo Fin de Máster Máster de Diseño Avanzado en Ingeniería Mecánica Produccion masiva de microburbujas monodispersas para aplicaciones reales evitando la microfluídica Autor: Enrique Jesús Sánchez Quintero Tutor: José Manuel Gordillo Arias de Saavedra Departamento Mecánica de Fluidos e Ingeniería Aeroespacial Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2017 Trabajo Fin de Máster Máster de Diseño Avanzado en Ingeniería Mecánica Produccion masiva de microburbujas monodispersas para aplicaciones reales evitando la microfluídica Autor: Enrique Jesús Sánchez Quintero Tutor: José Manuel Gordillo Arias de Saavedra Catedrático de Universidad Departamento Mecánica de Fluidos e Ingeniería Aeroespacial Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2017 Trabajo Fin de Máster: Produccion masiva de microburbujas monodispersas para aplicaciones reales evitando la microfluídica Autor: Enrique Jesús Sánchez Quintero Tutor: José Manuel Gordillo Arias de Saavedra El tribunal nombrado para juzgar el trabajo arriba indicado, compuesto por los siguientes profesores: Presidente: Vocal/es: Secretario: acuerdan otorgarle la calificación de: El Secretario del Tribunal Fecha: Agradecimientos Es te trabajo no hubiera sido posible sin la colaboración de mi tutor José Manuel Gordillo, a quien debo la gratitud de haber podido disfrutar con la realización de este proyecto. La motivación trasmitida durante todo el trabajo experimental y su rigurosidad en cada razonamiento me han hecho crecer tanto académica como personalmente. De la misma forma puedo estar agradecido a Álvaro Evangelio, quien comenzó este peculiar proyecto y cuyo trabajo ha sido capital para lograr los resultados aquí mostrados. También quiero agradecer la colaboración de compañeros como Manuel González en la fabriación y diseño de todos los dispositivos empleados en los experimentos y a Francisco del Campo por su ayuda con el montaje experimental. Como suele decirse, "y por último aunque no por ello menos importante", a mi pareja y mi familia, por todo su apoyo mostrado durante estos meses de trabajo. !Muchas gracias a todos! Enrique Jesús Sánchez Quintero Sevilla, 2017 I Resumen La producción de microburbujas monodispersas hoy en día ha alcanzado un importante desarrollo. Por su parte, las aplicaciones a nivel científico, industrial y médico que demandan cada vez más la generación de burbujas de tamaños submilimétricos han contribuido a motivar un gran incremento en el conocimiento de este campo de la Mecánica de Fluidos. Tanto es así que se han desarrollado multitud de tecnologías basadas en dispositivos microfluídicos con el fin de conseguir producir burbujas monodisperas con el menor diámetro posible y a las mayores frecuencias de burbujeo posibles. Si bien cada tipo de dispositivo difiere ligeramente de los demás y presenta sus ventajas e inconvenientes frente a estos, bien es cierto que todos tienen en común las micrométricas dimensiones características de su geometría. Este hecho implica que desarrollar y operar un dispositivo para la generación de microburbujas monodisperas no sea una labor sencilla; además, el empleo de dispositivos de esta escala limita mucho las frecuencias máximas de producción de burbujas. Sin embargo, el proceso de formación de una burbuja no va necesariamente ligado al empleo de este tipo de dispositivos microfluídicos. El estado del arte sobre el conocimiento en lo que al proceso de formación se refiere identifica, entre otras variables, al gradiente favorable de presión local en el punto de formación de la burbuja como uno de los actores principales en el proceso, determinando el diámetro y frecuencia de producción. Así, dado que existen multitud de geometrías donde un líquido puede verse sometido a fuertes gradientes favorables de presión, se ha decidido explorar otras soluciones tecnológicas que, empleando los mismos principios físicos que gobiernan la formación de burbujas en los dispositivos microfluídicos, permitan producir estas en geometrías completamente diferentes, con el fin de extender los campos de aplicación de la producción monodispersa de microburbujas y aumentar las frecuencias actuales de producción. III XÍndice Apéndice A Introducción a la Aerodinámica potencial 45 A.1 Introducción 45 A.2 Ecuaciones generales 46 A.3 Ecuaciones simplificadas 47 A.4 Ecuación de la vorticidad y flujos potenciales 51 A.5 Resolución de la ecuación de Laplace mediante superposición de soluciones elementales 52 A.5.1 Función potencial correspondiente a una corriente uniforme, φ152 A.5.2 Función potencial correspondiente a un manantial o un sumidero, φ253 A.5.3 Función potencial correspondiente a un torbellino, φ354 A.5.4 Función potencial correspondiente a un doblete, φ454 A.5.5 Potencial creado por la superposición de soluciones elementales 55 A.6 Fuerzas sobre objetos bidimensionales según la teoría potencial 56 A.6.1 Fuerzas sobre un cilindro de radio r=ainmerso en un flujo potencial 56 A.6.2 Fuerzas sobre objetos bidimensionales de geometría arbitraria inmersos en un flujo potencial 57 A.7 Generación de circulación y condición de Kutta-Joukowski 58 A.8 Resolución numérica del flujo alrededor de un perfil aerodinámico: método de los elementos de contorno 61 Apéndice B Aerodinámica de alas de envergadura finita. Vortex-Lattice 67 B.1 Ecuación del potencial y condiciones de contorno 67 B.2 Resolución del problema sustentador: método de Vortex Lattice 70 B.2.1 Campo de velocidades. Analogía de Biot y Savart 70 B.2.2 Aspectos relevantes de la implementación numérica 73 B.2.3 Circulación y ángulo de ataque inducido 74 Índice de Figuras 77 Índice de Tablas 81 Bibliografía 83 Índice alfabético 85 Glosario 85 1 Introducción Pr oducir burbujas de tamaño micrométrico tiene, hoy día, numerosas aplicaciones que van más allá de aquellas destinadas a procesos a escala de laboratorio. De hecho, el gran desarrollo y la diversidad de tecnologías desarrolladas en los últimos años han dado lugar a recientes revisiones del estado del arte que tratan con gran detalle tanto los fundamentos que sustentan la producción de microburbujas monodispersas como sus aplicaciones (véase [ 1 ]). Tanto es así que el empleo de microburbujas tiene cabida en procesos de índole tan variada como el tratamiento de aguas, la industria alimentaria, y todo tipo de procesos médicos y farmacológicos, por citar algunos ejemplos. Así, para la obtención de imágenes por ultrasonidos, el uso de microburbujas como angentes de contraste ha demostrado arrojar unos excelentes resultados [ 2 , 3 , 4 ]. Por otro lado, las elevadas necesidades de aireación y el gran porcentaje que esta ocupa en el coste de operación de los biorreactores [ 5 , 6 ] hacen que el adecuado control del tamaño y la frecuencia de producción de las burbujas tenga un fuerte impacto en la eficiencia del proceso. Cuando se requiere el uso de microburbujas en aplicaciones como las mencionadas anteriormente, se dispone de tres variables que se desean controlar: el diámetro medio de las burbujas, db , la frecuencia de producción de estas, fb y el índice de polidispersión, PDI =∑n i(Di−hDi)2/n (siendo Di el diámetro de cada burbuja, hDi el diámetro medio y n el número de partículas), ya que para considerar que la producción de microburbujas es monodispera estas tienen que tener un PDI por debajo del 5% [ 1 ]. Así, en el caso de la oxigenación de biorreactores, típicamente será necesario satisfacer la demanda de oxígeno de los microorganismos presentes, OUR (Oxygen Uptake Rate, de sus siglas en inglés), por lo que el ratio de transferencia de oxígeno (Oxygen Transfer Rate-OTR) puede ser el paso limitante en todo el proceso [ 5 ]. El parámetro que controla el OTR es, dado un gradiente de concentración, el coeficiente de transferencia de masa, kLa , siendo kL el cociente entre la constante de difusión para una diferencia de concentración y área determinadas, y se encuentra direcamente afectado por la frecuencia y el diámetro de las burbujas. En efecto, el area específica de interfase, i.e. el area por unidad de volumen de las burbujas, depende de forma inversa del diámetro de las burbujas ( a=6φ/db , con φ la fracción de gas en el medio). Por lo tanto, a menor diámetro de las burbujas mayor es el coeficiente de transferencia de transferencia de masa, lo que se traduce en última instancia en una reducción del caudal de aire requerido para satisfacer el OUR del cultivo, con el consiguiente ahorro que esta implica. Es por ello que el empleo de difusores de burbuja fina [ 7 , 6 ], con diámetros 1-3 mm, o incluso de microburbujas [ 8 , 9 , 10 ], con diámetros 10-500 µ m, resulta esencial cuando se trata de aumentar la eficiencia del proceso de aireación; no obstante, también se ha reportado que la presencia de surfactantes y antiespumantes en el medio puede reducir críticamente el valor de dicho coeficiente en los casos en los que se favorece la coalescencia (véase [ 11 ]), mientras que la presencia de una determinada concentración de sal (equivalente al lodo presente en las aguas no tratadas) puede inhibir precisamente esta coalescencia [7]. Como puede observarse, las exigencias de las industrias actuales, que requieren cada vez mayores frecuencias de producción y menores diámetros de las burbujas, conlleva que se hayan desarrollado diferentes tecnologías para poder conseguir una población monodispersa donde se pueda controlar tanto db como fb . Todas estas tecnologías emergentes emplean dispositivos microfluídicos de tamaño milimétrico y submilimétrico que, aunque puedan parecer muy similares entre sí, se fundamentan en principios físicos diferentes que conviene comprender [ 1 ]. De este modo, el capítulo se estructura de la siguiente forma: en primer lugar, se realiza una descripción de las ecuaciones que gobiernan la dinámica de una burbuja que se produce en el seno de un líquido, después, se enumeran las diversas tecnologías que se han desarrollado y que se utilizan para generar una población monodispersa de microburbujas en los diferentes regímenes. En algunas de esta aplicaciones, el papel del gradiente de presión existente en el líquido exterior juega un papel fundamental 1 2Capítulo 1. Introducción que se describe con detalle en la Sección 1.3. Finalmente, considerando las pequeñas geometrías empleadas en las tecnologías actuales y el rol determinante del gradiente de presión, se expondrán las ideas que han llevado al desarrollo del dispositivo que en este trabajo se estudia. 1.1 Fundamentos de la generación de burbujas Lejos de lo que pudiera intuirse, la generación de burbujas es un proceso con importantes diferencias respecto al procceso de producción de gotas. En general, para generar gotas de radio rd a una frecuencia fd∼U/rd , basta con inyectar el líquido a través de un tubo de radio rt∼rd a una velocidad U&Uc , con Uc= [σ/(ρg)]1/2 la velocidad capilar1[1]. Sin embargo, en el caso de la formación de burbujas para estas condiciones, el comportamiento es diferente. Por un lado, en el caso de generación cuasiestática en el que UUc , el radio de la burbuja viene dado por el conocido radio de Fritz, que resulta del balance de esfuerzos de tensión superficial con la fuerza de flotación, rF/rt∼[3/(2Bo)]1/3, con Bo =ρgr2 t/σel número de Bond, que mide la importancia relativa de los esfuerzos de tensión superficial frente a los esfuerzos de volumen (gravedad/flotación). Sin embargo, si se desea reducir el diámetro de las burbujas o incrementar la frecuencia de producción aumentando la velocidad de inyección del gas a valores por encima de Uc , lejos de obtener un chorro de radio comparable al del injector, se obtienen (por encima de cierta velocidad) burbujas de volumen Vb∝Qg/g1/26/5>4/3πr3 F [ 1 ], y lo que es más, si se sigue aumentando la velocidad, las burbujas cercanas entre sí pueden coalescer, con lo que el diámetro final obtenido es mucho mayor [?]. Para explicar estas diferencias entre ambos procesos, si se desprecian los efectos dinámicos del gas (ya que ρg/ρ1yµg/µ1), y se considera que la burbuja es prácticamente esférica, con una presión uniforme en su interior, la dinámica de la burbuja puede describirse a través de la ecuación de Rayleigh-Plesset. ρRb¨ Rb+3 2˙ R2 b=∆pexit −2σ Rb−4µ˙ Rb Rb (1.1) Tal y como se puede observar en [ 1 ], para que el crecimiento de la burbuja y su posterior colapso tengan lugar, requiere que el término ∆p−2σ/Rb cambie de positivo, valor que adquire en los primeros instantes mientras la burbuja se infla, a negativo, momento en el que las velocidades negativas cerca del inyector la harán colapsar. La Ecuación (1.1) junto con la ecuación de continuidad Qg=dVb dt≃4πR2 b˙ Rb(1.2) y la relación entre la presión de inyección y el caudal en la línea de gas2 p0−pexit =ρgKRegQg/πr2 t2(1.3) proporcionan el conjunto de ecuaciones necesario para describir sucintamente las tecnologías de la Sección 1.2 y, con ellas, el diámetro y frecuencia de producción de las burbujas finalmente obtenidas. Sintetizando una vez más las ideas expuestas en [ 1 ], el proceso de formación de una burbuja en una piscina en reposo puede esquematizarse en las siguientes etapas: • En primer lugar, para satisfacer Ecuación (1.3) , el volumen de la burbuja debe aumentar, lo que provoca velocidades radiales en el líquido hacia fuera de la burbuja. • La presión del gas en el interior de la burbuja es prácticamente uniforme, y debe adaptarse a la del líquido exterior en algún punto entre z=0 (véase la Figura 1.1) y z=2Rb , el cenit de la misma, supóngase en z∼Rb. • A medida que el radio de la burbuja aumenta, la punta de la burbuja se encuentra con una sobrepresión respecto al líquido exterior del orden de ∼ρgRb , lo que contribuye a que el gas acelere al líquido exterior, mientras que, en z=0 , es el líquido el que ejerce una sobrepresión sobre la base de la burbuja del orden de ∼ −ρgRb , lo que conlleva que el líquido induzca velocidades hacia el interior de la burbuja. 1Se remite al lector interesado en la rotura e inestabilidades de chorros líquidos a la excelente revisión [12] 2 Nótese que no puede presuponerse, a priori, la constancia del caudal, debido a las dos etapas bien diferenciadas que existen en el proceso de generación de burbujas. 1.2 Dispositivos para la generación de microburbujas 3 7 6 5 4 3 2 1 00 1 2 3 0 1 2 3 4 5 0 2 4 6 8 380 Pressure (Pa) r/rt z/rt 360 340 320 300 280 Figura 1.1 Contornos de presión en las fases iniciales del proceso de formación de una burbuja para un número de Bond, Bo =0.245. Cortesía de [1]. Por lo tanto, si nos ceñimos al caso no viscoso, la burbuja se desprenderá del inyector cuando la velocidad hacia adentro en z=0 , que surge del balance de Ecuación (1.1) , ∆p=pexit −2σ/Rb∼ρ˙ R2 b coincida con la velocidad hacia afuera impuesta por continuidad, ∼Qg/4πR2 b. Qg R2 b∼s∆p ρ∼pgRb⇒db∼Qg g1/22/5 (1.4) siendo, finalmente, la frecuencia de producción fb∼Qg d3 b∼ g3 Qg!1/5 (1.5) Así, aunque puede emplearse este sencillo método para producir burbujas simplemente inyectando gas en el seno de un líquido en reposo, la burbujas obtenidas poseen un diámetro significativamente mayor que el radio del inyector y con frecuencias que decrecen con el caudal, por lo que no resulta un método adecuado para satisfacer las demandas que las aplicaciones actuales requieren en lo que a diámetros y frecuencias se refiere[ 1 ]. Ello ha propiciado el desarrollo de tecnologías más sofisticadas basadas en dispositivos con geometrías microfluídicas como las que se detallan a continuación, capaces de producir en masa burbujas de diámetros micrométricos. 1.2 Dispositivos para la generación de microburbujas Una vez que se han descrito las ecuaciones necesarias para la compresión de los fundamentos de la generación de burbujas en el seno de un líquido, se está en disposición de enumerar y describir de forma sucinta las tecnologías más relevantes que, actualmente, se emplean para producirlas masivamente las mencionadas burbujas. Conviene recordar que las soluciones tecnológicas que aquí se presentan no son las únicas que permiten la producción de burbujas con tamaños submilimétricos; en efecto, en aplicaciones industriales como los reactores químicos, los esfuerzos de cortadura causados por el flujo turbulento generan burbujas micrométricas pero con una gran dispersión de tamaños [1]. De nuevo, al igual que se realizó en Sección 1.1, se respetará la estructura de [ 1 ] para presentar las tecnologías que a continuación se describen, no extendiénsose en exceso y pudiendo encontrar el lector una descripción más detallada tanto en el review como en las referencias citadas en él. Así pues, en la Figura 1.2, se describen de forma esquemática los dispositivos más relevantes para la producción de burbujas actualmente. En la figura se pueden distinguir dos tipos principales de tecnologías: aquellas emplean una corriente de líquido para provocar el colapso de la corriente gaseosa para la formación de burbujas (Figura 1.2a-d) y aquellas que en las que se emplean ultrasonidos (Figura 1.2e-f). A su vez, de entre las primeras, se puede distinguir aquellas en las que el líquido es inyectado en la misma dirección que la corriente gaseosa (Figura 1.2a-c) y aquellas en las que el líquido y el gas fluyen de forma perpendicular (Figura 1.2d). Finalmente, la diferencia principal entre los dispositivos de flow-focussing (Figura 1.2b-c), y los de coflow es que en el primero se hace 4Capítulo 1. Introducción U UQ g 2r t2r t p a p 0 Gas Liquid Simple co  ow a Cross  ow d Sonication device e Ultrasound bubbler f Axisymmetric  ow focusing b Planar  ow focusing c Q Q g Q g(t) Q g(t) D Q Q Q gQ Focusing channel exit channel w w exit U w Ultrasound wave Gas meniscus l f l Figura 1.2 Representación esquemática de los diferentes dispositivos que identiican la tecnología utilizada en la producción de burbujas monodispersas. Figura adaptada de [1]. pasar ambos fluidos a través de un estrechamiento, lo que acelera colapso y permite generar burbujas más pequeñas [ 1 ]. En esta sección, nos centraremos sólo en los 4 primeros, debido a la analogía que los mismos presentan con la solución tecnológica de la que trata este estudio. 1.2.1 Dispositivos de Coflow Un dispositivo de tipo coflow es aquél como el mostrado en la Figura 1.2a, donde la corriente gaseosa y la de líquido son inyectadas en la misma dirección. Las condiciones en las que usualmente se opera este tipo de dispositivos son aquellas en las que el número de Reynolds y Weber son tales que Re =ρUrt/µ1 y We =ρU2rt/σ . Bajo estas condiciones, las frecuencias de producción, fb y los diámetros equivalentes de las burbujas obtenidas, db=6Qg/(πfb)(1/3) , dependen sólo del ratio de velocidades gas-líquido, Ug/U , y del ratio rt/U . Además, en los casos bajo consideración donde We 1 , el líquido exterior impone la velocidad a la que la interfase es transportada, con lo que se previene la coalescencia ya que las burbujas son transportadas a la velocidad del líquido [ 13 ]. Precisamente, en una operación normal con este tipo de dispositivos, la velocidad del gas suele ser mayor que la del líquido, con lo que las burbujas son generadas cerca de la punta del inyector. Siguiendo las mismas ideas que el proceso basado en las fases de expansión y 1.2 Dispositivos para la generación de microburbujas 5 posterior colapso descrito en Sección 1.1, en [ 14 ] se desarrolla un modelo simple para la fase de colapso en el que puede comprobarse que, en el límite en el que Ug/U1, la frecuencia escala como fb≃0.15U/rt. Por otro lado, el caso contrario en el que Re 1 no ha sido muy reportado en la literatura [ 1 ]. En este caso, en lugar del número de Webber, el parámetro adimensional que gobierna el problema es el número capilar, Ca =µU/σ , además del ratio Ug/U como en el caso anterior y del ya descrito número de Bond, Bo . En este caso, tras el exhaustivo estudio numérico de [ 15 ], se tiene que el diámetro de las burbujas decrece cuando el ratio Ug/Utambién lo hace. 1.2.2 Dispositivos de Cross-Flow Los dispositivos de tipo Cross-Flow son aquellos como los mostrados en la Figura 1.2d, donde puede apreciarse que el gas es inyectado de forma perpendicular al líquido a través de una junta en T. Para tener buen control y repetibilidad en los tamaños de las burbujas, estos dispositivos operan a números de Reynolds muy bajos, por lo que su mayor campo de amplicación se encuentra dentro de la microfluídica. Dado que operan a Re 1 , se observan diferentes regímenes en función del número capilar, Ca , como son el squeezing, dripping ojetting. De entre ellos, el que parece más apropiado para el control del diámetro de las burbujas es el squeezing, lo que ocurre cuando Ca .O10−2 . Por otro lado, una de las principales desventajas del empleo de este tipo de dispositivo es la baja frecuencia de producción derivada de la operación a bajos números capilares. En efecto, si se pretende aumentar la velocidad del líquido para aumentar dicha frecuencia, también se producirá un aumento del numero capilar, lo que limita en general la frecuencia de producción de estos dispositivos a fb∼103Hz. Tanto los dispositivos con configuraciones de tipo Cross-Flow como los de Coflow permiten tener control de forma separada del diámetro de las burbujas y de su frecuencia, variando para ello Qg y Q . Sin embargo, a pesar de ser dispositivos microfluídicos de tamaños del orden de las centenas de micras(o mayores), su geometría limita el tamaño mínimo de las burbujas que se pueden obtener. Esta dificultad puede ser superada a través de otro método conocido como Flow-Focussing que utiliza una geometría un tanto diferente como se verá en la próxima subsección. 1.2.3 Dispositivos de Flow-Focussing La sencillez de los dispositivos descritos anteriormente posee un notable inconveniente, pues si se elimina la inyección de gas las líneas de líquido son casi paralelas, con lo que los gradientes de presión existentes son prácticamente despreciables [ 1 ]. La demanda de diámetros cada vez menores de las burbujas hace que sea tecnológicamente inviable reducir tanto como se desee el diámetro del inyector, por lo que otras soluciones han emergido para paliar las limitaciones de las configuraciones de coflow y crossflow; una de ellas, es el empleo de dispositivos de Flow-Focussing. En geometrías de flow-focussing como las mostradas en la Figura 1.2b-c, la presencia del estrechamiento provoca (además de un coflujo de líquido y gas que previene la aparición de coalescencia) un gradiente longitudinal de presión, cuyo efecto principal es el de reducir el diámetro final de las burbujas obtenidas. Los primeros en reportar este fenómeno fueron [ 16 ], obteniendo burbujas de tamaños db∼O(10µm)<D , controlando simplemente el ratio Qg/Q y el caudal de líquido exterior Q=πD2U/4 , con D el diámetro del estrechamiento. En efecto, el papel del gradiente axial de presión es análogo al jugado por la gravedad en el caso de la Ecuación (1.4) , esto es, tomando ∆pexit =∇pRb , y teniendo en cuenta que el gradiente de presión puede ser obtenido de forma estimada de las ecuaciones de Navier-Stokes como ρDu Dt =−∇p⇒∇p∼ρu·∇u∼ρU2/D(1.6) se tiene, sustituyendo ∆pen (1.4), db D∼Qg Q2/5 (1.7) siendo pues la frecuencia de burbujeo fb=6Qg πd3 b ∝U DQg Q−1/5 (1.8) por lo que este dispositivo también permite controlar de forma separada tamaño y frecuencia de producción simplemente a través de la variación de los ratios Qg/Q y de la velocidad del líquido U . Uno de los principales 6Capítulo 1. Introducción escollos en el uso de este tipo de dispositivos se encuentra en la dificultad que supone conseguir un buen alineamiento entre los diferentes canales, lo que introdujo la configuración de la Figura 1.2c, que permite simplificar esta labor. 1.2.4 Otros dispositivos y breve reflexión Existe, además de los descritos en los apartados anteriores, otros dispositivos para generar burubujas de tamaño micrométrico de forma monodispersa, como son aquellos basados en la aplicación de ondas acústicas, tal y como se muestra en la Figura 1.2e-f, y lo que es más, nada impediría (en principio) combinar estos dispositivos con los anteriores, si se pretende reducir aún más el diámetro de las burbujas [ 1 ]. A través de los ejemplos más significativos de dispositivos que aquí se han presentado, el lector puede estimar que el tamaño característico de estos dispositivos se encontrará en un rango comprendido entre ∼O(10µm) hasta ∼O(10mm) , lo que a priori puede suponer una dificultad tanto de fabricación como de operación; además, a la dificultad de manejo de los instrumentos microfluídicos habría que sumar la limitación existente en la frecuencia máxima de producción, que estaría limitada por el tamaño de los dispositivos y cuya solución se encontraría en la fabriación de un número mayor de estos. Cabe entonces preguntarse si los mecanismos físicos ya descritos y que subyacen bajo los diferentes procesos de producción de microburbujas aquí presentados son exclusivos de los dispositivos microfluídicos. En efecto, si bien a lo largo del capítulo se ha podido comprobar que existe una clara dependencia con el diámetro del inyector del gas (lo que apoya la necesidad del empleo de la microfluídica), también es cierto que las condiciones de contorno que el líquido exterior impone sobre la corriente gaseosa, que son que la entrefase se convecte aproximadamente a la velocidad del líquido (coflow y cross-flow) y el posible gradiente longitudinal de presión provocado por un estrechamiento (flow-focussing), no tendrían (en principio) porqué ir ligadas taxativamente al empleo de geometrías micrométricas. Así pues, en la siguiente sección de este capítulo se discute con más detalle el papel que el gradiente de presión juega en la formación de burbujas, prestando especial atención para ello en una técnica conocida como Confined-Selective-Withdrawal que presenta características comunes tanto a las configuraciones de coflow y flow-focussing. Se pretende así reforzar la idea de la importancia que juega el gradiente de presión en la formación de burbujas para, finalmente, plantear si cabría la posibilidad de diseñar un dispositivo capaz de generar microburbujas monodispersas de forma masiva cuya dimensión geométrica característica fuera mucho mayor que las de los dispositivos empleados actualmente. 1.3 Influencia del gradiente de presión A lo largo de la Sección 1.2, se ha podido comprobar como el papel del líquido exterior en la generación de microburbujas puede ir desde simplemente imponer la velocidad de la entrefase previniendo la coalescencia hasta ejercer un efecto análogo al de la gravedad a través del gradiente axial de presión. Las conclusiones extraídas en la Subsección 1.2.4 motivan indagar un poco más en el papel que dicho gradiente longitudinal de presión puede jugar en el proceso de formación de burbujas, por lo que se van a presentar de forma sucinta los resultados obtenidos [ 17 ] donde se obtuvieron burbujas de tamaño micrométrico empleando la técnica conocida como Confined Selective Withdrawal. En la Figura 1.3 se describe esquemáticamente el funcionamiento de un dispositivo basado en esta técnica. Consiste en cámara presurizada a p0 con un líquido, de densidad ρ y de viscosidad µ , que se suministra desde el exterior y que descarga a pa<p0 a través de un tubo capilar (en este caso de sección cuadrada) de dimensiones características L∼1mm . El capilar se encuentra perfectamente alineado con una aguja donde se inyecta el gas con caudal Qg y que se encuentra a una distancia H del capilar. Como puede observarse en la Figura 1.4, para velocidades del líquido en el capilar inferiores a un determinado valor, U<U∗ , el menisco emite burbujas de forma pulsante y de diferentes diámetros [ 17 ]. Sin embargo, por encima de esta velocidad, se puede observar que el menisco de aire se mantiene estable y estacionario en una regíon x≤xs (véase la Figura 1.5), lo que propicia la aparición de un régimen de producción de burbujas monodispersas; la condición de equilibrio de la que resulta xs será discutida más adelante. También, en la Figura 1.5 se muestra como, para valores de U>U∗ , aguas abajo de la zona estacionaria ( x≤xs ), se produce un cilindro de gas de diámetro dg y longitud ` que, finalmente, culmina con la producción de una nueva burbuja; una vez la burbuja es emitida y transportada a la velocidad del líquido, se inicia nuevamente el proceso de formación de otra burbuja a través del mecanismo descrito. Aunque en [ 17 ] se analiza tanto el caso viscoso como el no viscoso, en esta sección nos centraremos únicamente en los casos Re 1 , ya que será el caso de interés para el dispositivo que aquí se describe. Antes de comenzar el análisis, conviene notar que la presión del gas, pg , que además se considera constante a lo 1.3 Influencia del gradiente de presión 7 H L 1 Gas injection needle1 Extraction tube (square duct) 2 2 Bubble 3 3 4 5 5Pressurized chamber6 6 Liquid inflow7 7 4End of the extraction line Air meniscus Figura 1.3 Esquema del dispositivo empleado en [ 17 ] para la producción de microburbujas con la técnica de Confined Selective Withdrawal. Imagen adaptada de [17]. Fully unsteady gas ligament Steady portion of the gas ligament Unsteady portion of the gas ligament where bubbles are formed Needle tip (a) (b) Figura 1.4 Diferentes regímenes de burbujeo en función de la velocidad U a la que el líquido circula a través del capilar. En (a) se visualiza el caso en que U<U∗ , por lo que no se tiene un régimen de burbujeo constante, mientras que en (b), U>U∗ , se consigue alcanzar un régimen de producción de burbujas monodispersas. Imágenes adaptadas de [17].. largo todo el menisco según la hipótesis ya comentada en la Sección 1.1, debe de encontrarse en equilibrio con la del líquido a la salida de la aguja, donde la velocidad del líquido es casi nula, por lo que puede suponerse que pg≃p0 . Así, para obtener un menisco de aire estable en la región x≤xs , la diferencia de presiones p0−p(x) debe estar en equilibrio con la presión capilar. En efecto, considerando que el hilo gaseoso tenga un diámetro dg , llamando U0 a la velocidad en el centro del tubo capilar, y teniendo en cuenta el cumplimiento de la ecuación de continuidad Qg=πd2 g 4U0⇒dg∼Qg U01/2 ,(1.9) 8Capítulo 1. Introducción se tiene que p0−p(xs) = 2σ dg≃2σ Qg/U01/2(1.10) lo que proporciona la condición de equilibrio para que exista un menisco de aire estacionario en x≤xs . Tanto en (1.9) como en (1.10) , la expresión de U0 vendrá determinada por el tipo de flujo que exista en el capilar, lo que para el caso en el que Re 1 y teniendo en cuenta los resultados de las simulaciones 3 realizadas en [ 17 ] constituye un perfil uniforme de velocidad; efectivamente, este hecho queda representado por el valor del coeficiente de presión adimensional ξ= (p0−p(x))/1/2ρU2≃1 en x=0 . Así, dado que p0−p(xs)∼1/2ρU2 para Re 1 y teniendo en cuenta que el proceso de formación de burbujas será poisble si p0−p(x)∼1/2ρU2 0>2σ dg , son condiciones necesarias para la producción de burbujas monodispersas que β=ρU2L 4σq1/2&1 conq=Qg UL2(1.11) donde β en la Ecuación (1.11) constituye un parámetro similar al número de Weber, We =ρU2L2/σ ; la validez de (1.11) queda mostrada en [17] Steady Unsteady (a) (b) Figura 1.5 Secuencia de producción de burbujas a partir de un menisco de aire estable para la región x≤xs . Como se aprecia en la figura, aguas abajo de xs , se emite un cilindro de gas de diámetro ∼dg que se extiende una longitud `. Imagen tomada de [17]. Una vez se han descrito las condiciones neesarias para la producción monodispersa de burbujas y se conocen las condiciones de equilibrio estacionario del menisco de aire, se está en condición de describir el resto del proceso de formación, así como las expresiones de los diámetros y frecuencias de formación de las burbujas obtenidas. Una vez que las burbujas se desprenden del ligamento de gas en Lx ≈Lxs+` , la diferencia de presión ∆p=p0−p(x)−2σ/dg>0 , induce velocidades radiales sobre el ligamento de gas que provoca el inicio del proceso de formación de una nueva burbuja. Dado que durante los instantes posteriores a la eyección de una nueva burbuja el diámetro del cilindro apenas varía, es posible escribir ∆p=p0−p(x)−2σ/dg=p0−p(xs)−2σ/dg+p(xs)−p(x)≈ ≈−dp dx(xs)(x−xs)≈d(p0−p) dx(xs)` L (1.12) 3 En [ 17 ] se realiza una simulación del dominio considerando este como axilsimétrico y sin simular la fase gaseosa, lo que suele realizarse generalmente cuando se desea estudiar la estabilidad de chorros, descomponiendo el campo de presiones como suma de la solución básica más una perturbación (véase como ejemplo [ 18 ]). En [ 17 ], sin embargo, se realiza la hipótesis, verificada a posteriori, de que la presencia de las burbujas no perturba el campo de presiones del líquido. 1.4 Analogía Aerodinámica 9 donde se ha tenido en cuenta la Ecuación (1.10) y se ha realizado un desarrollo en serie de Taylor de primer orden de la presión en torno al punto x=xs . De este modo, y tal y como se detalla en [ 17 ], la Ecuación (1.12) muestra como el proceso de formación de burbujas está gobernado por el gradiente de presión local en el punto x=xs , por lo que, como veníamos anunciando, el gradiente de presión local posee un papel absolutamente análogo al que tiene la gravedad en la formación de una burbuja en una piscina en reposo. Así, siguiendo el procedimiento que se siguió en la Sección 1.1, si se sustituye la Ecuación (1.12) , particularizada para el caso Re 1, en la ecuación de Rayleigh-Plesset (Ecuación (1.1)), se tendrá que ρRb¨ Rb∝d(p0−p) dx(xs)` L∝` 2LρU2Ps(1.13) siendo Ps=˙ ξ(xs). La ecuación anterior junto con la ecuación de continuidad para la línea de gas Qg=π 6d3 bfb(1.14) permiten obtener las frecuencias de producción, ρPs 2 U2 L`∝ρRb¨ Rb∼ρd2 bf2 b⇒fb∝UpPs/2 pLdbs` db (1.15) y, finalmente, tomando `∝db por analogía al caso de formación de burbujas en una piscina en reposo [ 17 ], la expresión final de la frecuencia y los diámetros de las burbujas. fb∝UpPs/2 pLdb ydb L∝ Qg UpPs/2L2!2/5 (1.16) El proceso seguido para la obtención de las ecuaciones (1.15) y (1.16) es completamente análogo al descrito en la Sección 1.1 con el gradiente de presión local ejerciendo el papel de la gravedad y con la diferencia que ρU2/Lρg , lo que atendiendo a las ecuaciones de arriba, se traduce en un aumento de la frecuencia y una disminución de los diámetros obtenidos. Este hecho invita a pensar que, si el proceso de formación descrito está controlado por el gradiente de presión local en xs , nada impide imaginar dispositivos generadores de burbujas donde los gradientes favorables de presión puedan obtenerse con geometrías muy diferentes a la descrita en [17]. 1.4 Analogía Aerodinámica Al final de la Sección 1.3, basándose en las ideas de [ 17 ], se dejó entrever la posibilidad de obtener gradientes favorables de presión similares a los encontrados en dispositivos de Confned Selective Withdrawal o Flow- Focussing pero con otras geometrías completamente diferentes. La motivación en la exploración de otras geometrías es doble: dispositivos con tamaños alejados de la escala de la microfluídica permitirían no sólo una fabricación y operación más sencilla, sino también aumentar notablemente el ritmo de producción de las burbujas. Así, un ejemplo cotidiano de geometrías donde se producen fuertes gradientes favorables de presión lo constituyen las superficies sustentadoras de los aviones. Dado que la única premisa que, a priori, debe cumplir una geometría alternativa a las empleadas actualmente sería conseguir gradientes de presión comparables a los creados en los dispositivos microfluídicos, ¿qué impide pensar que se pueda emplear un ala para producir masivamente microburbujas monodispersas? Conviene antes de seguir, no obstante, esbozar algunas ideas básicas en aerodinámica que permitirán comprender mejor la posibilidad de emplear un ala como dispositivo generador de burbujas. Se ha decidido incluir en el Apéndice A un resumen de [ 19 ] donde se desarrollan los conceptos más importantes de la Aerodinámica potencial, de forma que el lector interesado pueda consultarlo para profundizar un poco más en las ideas y definiciones expuestas a continuación. La Aerodinámica es la parte de la Mecánica de Fluidos que se encarga del estudio del movimiento de gases (generalmente aire) alrededor de un cuerpo. Teniendo en cuenta las propiedades del aire (densidad ρg∼O(1kg/m3) y viscosidad µ∼O(10−5Pa ·s) a T=20◦ ) y las velocidades características de un cuerpo que se mueve en él (tómese como ejemplo un coche circulando a V=30m/s o el ala de un avión a V=100m/s ) se tiene que los números de Reynolds característicos del flujo de aire alrededor de un objeto con dimensiones características L=1m serán Re =ρVL/µ∼O106 . Bajo estas condiciones los efectos viscosos pueden ser despreciados en la mayor parte del dominio del problema 16 Capítulo 2. Ala acuática bidimensional Tabla 2.1 Dimensiones relevantes del túnel hidrodinámico. Dimensión Valor [cm] Ancho de la sección 30 Altura de la sección 30 distancia entre soportes 20 Figura 2.1 Perspectiva del túnel hidrodinámico donde los experimentos han sido realizados. en [ 17 ], esto es db∼O(200 µm) , y se estima que las burbujas son convectadas aguas abajo del perfil a una velocidad 2v∼U∞ , se tendrá que una burbuja recorre una distancia del orden de su diámetro en un tiempo tb∼db/U∞∼O(200 µs) , por lo que serán necesarios del orden de 5000 fotogramas por segundo (en adelante fps): resulta evidente que es necesario utilizar una cámara de alta velocidad. Tunel hidrodinámico En el laboratorio de Mecánica de Fluidos de la Escuela Superior de Ingeniería de Sevilla se dispone de un túnel hidrodinámico como el mostrado en la Figura 2.1. Como puede observarse en la figura, el túnel consiste en un conducto de sección rectangular a través del cual se hace circular un líquido (en este caso agua) utilizando para ello una bomba. Una vez que la bomba ha hecho ascender al líquido hasta la altura adecuada, se le hace pasar a través de una cámara de contracción cuyo fin es uniformizar el flujo incidente en la cámara de ensayos (zona transparente del túnel). Pese a no apreciarse en la imagen, el túnel se encuentra abierto y en contacto con la atmósfera en su parte superior, mientras que la parte inferior del mismo también esta delimitada por cristales de igual transparencia que las paredes laterales. Las dimensiones del túnel pueden consultarse en la Tabla 2.1 El caudal que circula a través del túnel hidrodinámico depende del producto de la altura total del agua (esto es, la fracción de la sección del túnel que queda mojada) y de la velocidad media del agua en la sección. Dado que el túnel se encuentra abierto en su parte superior, se ha comprobado experimentalmente que si la velocidad del túnel excede el valor de U∞∼0.75m/s se forman ondas superficiales. Por lo tanto, se tomará como restricción al diseño que la velocidad de la corriente incidente al perfil no puede superar este valor, es decir, U∞≤0.75m/s(2.1) Finalmente, cabe destacar que el control de la velocidad del túnel se realiza a través de un regulador de potencia, por lo que la medida de la velocidad para cada potencia debe determinarse experimentalmente. 2 En [ 17 ], las burbujas son convectadas a la velocidad de la corriente, siendo en este caso U la velocidad media en el capilar para Re 1. 2.1 Diseño y fabricación 17 Figura 2.2 Cámara de alta velocidad Phantom®v710. Reguladores de presión La formación de burbujas requerirá que se inyecte aire en ciertos puntos de la superficie del ala, el cual debe proceder desde dentro de la misma. Para poder suministrar el aire se dispone en el laboratorio de una red capaz de suministrar más de 6 bares de presión, lo cual es más que suficiente para el propósito que nos ocupa. Para regular de forma adecuada la presión se dispone de manorreductores cuya precisión es de 1 mbar. Cámara de alta velocidad Para la grabación de las imágenes se dispone de una cámara Phantom v710 como la que se muestra en la Figura 2.2, capaz de grabar a más de medio millón de fps. Dado que, según se comentó al inicio de esta sección, serán necesarios en torno a 5000 fps, esta cámara es capaz de acometer dicha tarea sin problema alguno. Adicionalmente, se disponen de dos objetivos Canon ® diferentes con zoom 1-4 y 2.5-10, respectivamente. Finalmente, cabe destacar que la distancia focal de la cámara es de aproximadamente 11 cm. 2.1.2 Diseño del ala Una vez que se consideran los equipos disponibles y las restricciones que pueda imponer cada uno, es momento de pasar a especificar el diseño final del ala. Dado que la sección transversal del ala consiste en un perfil aerodinámico, se ha decidido que el proceso de fabricación sea mediante control numérico, empleando un polímero de ABS (Acrilonitrilo Butadieno Estireno, de sus siglas en inglés) como materia prima. La fabriación por control numérico aúna las ventajas de la impresión 3D y la fabricación tradicional, ya que las paredes de la pieza obtenidas son siempre macizas 3 lo que aportará la robustez necesaria para soportar la experimentación del ala en el agua del túnel hidrodinámico. Para la sección transversal del ala se disponen de dos opciones principales en lo que a perfiles aerodinámicos se refiere: o bien se opta por diseñar un perfil personalizado, o bien se emplea alguno de los perfiles normalizados ya existentes, por ejemplo los de la serie NACA. Se ha elegido para este primer prototipo el conocido perfil simétrico NACA 0012, que posee un espesor máximo relativo a la cuerda t/c=0.12 . El motivo de empleo de un perfil de estas características se basa principalmente en que el NACA 0012 es uno de los perfiles aerodinámicos más estudiados, lo que unido a la sencillez en su diseño lo hacen idóneo para una primera prueba de concepto. No obstante, el empleo de un perfil con curvatura no supondría una pérdida de generalidad en los resultados obtenidos más adelante. Decidido el perfil aerodinámico que se empleará como sección transversal del ala, tan sólo es necesario definir la cuerda, c , y la envergadura, b , para tener completamente definida la geometría exterior del ala. Para la cuerda es necesario tener en cuenta los puntos de sujeción que tendrá el ala dentro del túnel hidrodinámico. En principio, podría optarse por disponer de 2 únicos puntos de sujeción (uno a cada lado del ala) que fueran acoplados en torno a una distancia c/2 del borde de ataque. Sin embargo, dado que en los experimentos se pretende variar el ángulo de ataque de la corriente inicidente, el diseño de un sistema capaz de variar la inclinación del ala con sólo 2 puntos de sujeción es bastante complejo. Por ello, resulta más sencillo disponer de 4 puntos de sujeción, de forma que se sitúen 2 de ellos cerca del borde de ataque y otros dos cerca del borde de salida, permaneciendo los primeros fijos durante la operación experimental, lo que permite variar la 3 Mediante impresión 3D también es posible emplear un 100% de entramado, con lo que se obtendrían resultados parecidos si se emplea el mismo material; si bien, el empleo de la técnica de impresión 3D con entramados tan elevados suele llevar aparejado un alto coste. 18 Capítulo 2. Ala acuática bidimensional inclinación del ala y con ello el ángulo de ataque, α , con el que inicide la corriente. El valor que toma la cuerda debe ser tal que permita que el ala se sitúe entre 2 puntos de soporte consecutivos del túnel (véase la Tabla 2.1), por lo que se ha decidido darle un valor c=30cm . En este momento, elegida la magnitud de la cuerda, se está en condiciones de poder calcular el número de Reynolds global en el ala, Rec=ρU∞c µ∼O1051.(2.2) donde se ha tomado un valor orientativo U∞∼0.5m/s . Por lo tanto, la magnitud de la cuerda escogida permite sin mayor complicación repdroducir las condiciones de Rec1 tal y como se pretendía. El cálculo del valor adecuado que se debe dar a la envergadura, por otro lado, sigue un razonamiento diferente. Modelar experimentalmente un problema bidimensional no es una tarea fácil en la mayoría de los campos de la física, ya que se deben asegurar el cumplimiento de las hipótesis que con el fin de simplificar el problema base se hayan realizado. Como se comentó en la Sección 1.4, para el caso del movimiento de un sólido fuselado en el seno de un fluido incompresible a altos números de Reynolds siendo el flujo inicialmente irrotacional, se obtiene que el problema se reduce al descrito por las ecuaciones (1.20) . La linealidad de este problema, admite que pueda ser resuelto mediante superposición de soluciones elementales de la laplaciana, como son términos fuente (y/o sumideros), torbellinos y dobletes; un ejemplo ilustrativo donde se emplea la superposición de soluciones elementales de tipo doblete lo consitiuye el Método de Green de Elementos de contorno descrito en el Apéndice A.Por otro lado, cuando existen contornos rígidos como una pared, la teoría de campos nos brinda el conocido método de las imágenes, consistente en aplicar una simetría en la distribución de soluciones elementales respecto a la pared en cuestión de forma que se cumplan las condiciones de contorno en la misma(en este caso, ∇φ·n=0 ). Este método permite modelar el efecto que la pared ejerce sobre el campo potencial solución del problema respecto al campo original. Este mismo razonamiento puede ser empleado de forma inversa: para modelar un perfil aerodinámico bidimensional, se tendría que disponer de un ala infinitamente larga, de modo que los torbellinos de la punta no indujeran velocidades transversales a la corriente incidente 4 y por tanto perpendiculares al plano del problema 2D, rompiendo así la hipótesis de bidimensionalidad. Gracias al método de las imágenes, es posible modelar este comportamiento si se empotra el ala en las paredes del túnel, para lo cual la envergadura del ala tendría que coincidir con el ancho de la sección. Sucede ,sin embargo, que las paredes del túnel son de cristal, lo que imposibilita la opción de atornillar o taladrar en ellas la sujeción del ala, de modo que se empleará una envergadura ligeramente inferior que permita un sistema alternativo de sujeción; nótese que, al no poder aplicar estrictamente el método de las imágenes, se está realizando pues la hipótesis de comportamiento bidimensional, que tendrá que ser corroborada a posteriori. El valor empleado de la envergadura es finalmente b=26cm. Hasta este punto se tiene completamente definida la geometría exterior del ala, pero aún no posee las condiciones para poder generar burbujas. Como ya se comentó previamente, el aire debe inyectarse a través de determinados orificios de la superficie del ala, por lo que el gas debe provenir directamente desde dentro de la misma. Para ello, existen dos posiblidades. Una sería llevar el aire desde un depósito situado en el exterior del túnel a través de conductos (como son los tubos de tygon o tubos de peek) hasta cada orificio, lo que sería directo pero complicaría el diseño, al tiempo que requeriría un sellado muy meticuloso de la zona de intersección de los tubos con la superficie exterior del ala. La otra, sería convertir el propio ala en un depósito, para lo cual simplemente habría que diseñar el dispositivo de forma que fuera hueco por dentro. Esta última opción es la que se ha elegido, pudiéndose observar en la Figura 2.3 un corte longitudinal en perpsectiva donde se aprecia la oquedad interior. Así, a través de un orificio roscado, un racor adaptará el interior del ala a la red de presión disponible. Para concluir el diseño del sistema de inyección de aire, es necesario indicar los puntos donde se situarán los orificios a través de los cuales el aire será introducido en el líquido. Recordando las conclusiones extraídas en el capítulo anterior, resulta evidente que los orificios deben situarse en las zonas donde el gradiente favorable de presión sea mayor, esto es, la región que, conteniendo al borde de ataque, se encuentra entre el punto de remanso y el pico de succión del perfil. No obstante, la posición de sendos puntos no es fija, sino que depende (para un perfil simétrico) del ángulo de ataque de la corriente incidente, por lo que parece conveniente distribuir los puntos de inyección en un entorno del borde de ataque y sobre una misma superficie del perfil (extradós o intradós), que será donde se grabe el proceso de formación de burbujas. De este modo, en cada punto y para cada α , se tendrán distintos valores del gradiente de presión adimensional, dCp/ds (véase la Figura 1.8a). Cabe apreciar además, que puede darse la circunstancia de que, para determinados valores de α , los puntos pasen de encontrarse en una zona de gradiente favorable de presión a una zona de 4Este aspecto pertenece a la teoría de alas de envergadura finita, cuya descripción detallada puede encontrarse en el Apéndice B 2.1 Diseño y fabricación 19 Figura 2.3 Perspectiva de la sección longitudinal del dispositivo generador de microburbujas. En la figura puede apreciarse el interior hueco del ala que se hallará presurizado para poder suministrar el caudal de gas requerido. Tabla 2.2 Posición de la abcisa, adimensionalizada con la cuerda del perfil, de los orificios de inyeccón de aire.. Puntos 1 2 3 4 5 6 7 8 9 x/c0 0.0017 0.003 0.005 0.0067 0.0083 0.01 0.02 0.03 gradiente adverso de presión, o lo que es lo mismo, aguas abajo del pico de succión. En la Tabla 2.2 se muestran las posiciones adimensionalizadas con la cuerda de la abcisa de los puntos de inyección tomando x según el eje longitudinal de simetría del perfil con origen en el borde de ataque. Finalmente, en lo que se refiere al diámetro y longtiud del los orificios, se debe tener en cuenta que resulta primordial asegurar la constancia del caudal durante el experimento, lo que no será trivial si se tienen en cuenta las perturbaciones que puedan existir tanto en la presión dentro del ala como en la presión en el otro extremo del orificio (al agua bordeando el perfil). Para minimizar las variaciones de caudal de gas causadas por las perturbaciones de presión y con idea de poder efectuar variaciones graduales de caudal mediante una variación mínima de presión de 1 mbar (precisión de los manorreductores), la solución que se propone es el empleo de tubos capilares que generen una gran pérdida de carga. En efecto, si se asume que el flujo en el interior del tubo capilar es tal que ReDD/Lt1 , 5 con ReD=4Qg/(πDν) el número de Reynolds basado en el diámetro del tubo D y Lt la longitud del mismo, el flujo de gas por el interior del tubo es casi-unidireccional y sigue la ley de Hägen-Poiseuille, Qg=πD4 128µg ∆p Lt (2.3) donde ∆p es el incremento de presón entre el interior y el exterior del ala. Si se llama ahora δp a la perturbación en el incremento de presión, se tiene que el incremento de presión total será ∆p=∆p0+δp , que provocará un flujo de air total igual a Qg=Qg0+δQg=πD4 128µgLt (∆p0+δp)⇒δQ=πD4 128µgLt δp=πD3 128µg D Lt δp(2.4) Por lo tanto, según la Ecuación (2.4) , se tiene que las perturbaciones en los incrementos de presión provacarán incrementos de caudal más pequeños cuanto menor sea el díametro D y el cociente D/Lt . De este modo, para que la inyección de aire a caudal constante sea lo más estable posible, se hace necesario reducir todo lo posible el diámetro del orifiio de salida. Así pues, considerando las limitaciones tecnológicas existentes para perforar un orificio del tamaño de, por ejemplo 200 µm , la mejor alternativa resulta en utilizar tubos capilares con el diámetro interior más pequeño posible compatible con un diámetro exterior más manejable, 5Esta hipótesis podrá ser comprobada una vez se analicen los resultados en la Sección 2.3 20 Capítulo 2. Ala acuática bidimensional (a) Perspectiva 1. (b) Perspectiva 2. Figura 2.4 Diferentes perspectivas del dispositivo generador de microburbujas bidimensional empleado en los experimentos. que permita insertarlos a través de los orificios del ala. Para el prototipo que se está describiendo, se han empleado tubos capilares de acero con diámetro exterior Dext =400 µmy diámetro interior D=160 µm. Para concluir esta sección, se muestran en la Figura 2.4 diferentes perspectivas del aspecto final del dispositivo generador de microburbujas descrito hasta este punto. Nótese la presencia en los extremos laterales del ala de un orificio pasante cerca del borde de ataque así como otro ciego cerca del borde de salida: estos orificios, que forman parte del sistema de variación del ángulo de ataque, serán descritos en detalle y por conveniencia dentro de la Subsección 2.1.3. 2.1.3 Diseño del montaje experimental Una vez ya se tiene completado el diseño y la fabriación del ala, tan sólo queda crear un montaje adecuado para poder realizar los experimentos. El primer objetivo de este montaje será un sistema capaz de no sólo posicionar el ala dentro del túnel hidrodinámico sino también de poder variar el ángulo de ataque, α . En la Figura 2.5 puede observarse que la estructura consta de los sigiuentes componentes: • Dos largeros longitudinales, los cuales sirven para anclar la estructura a los tornillos que sobresalen en la parte superior del túnel, de forma que la unión con este sea lo más robusta posible. • Dos largeros anteriores verticales, que ejercen la función de puntos de sujeción del ala cerca del borde de ataque. Estos largueros, que deben estar completamente en contacto con las paredes laterales del túnel, poseen orificios que permiten en paso del eje que recorre transversalmente al ala y que forma parte del sistema de variación del ángulo de ataque. • Dos largueros transversales, de los cuales el situado en la parte delantera (cerca del borde de ataque del ala) sirve para dar robustez a la estructura, mientras que el trasero (además de proveer estabilidad), incluye el sistema de variación del ángulo de ataque. El sistema de variación del ángulo de ataque, por lo tanto, funciona del siguiente modo: un eje metálico pasa a través del orificio pasante del ala cerca del borde de ataque (véase la Figura 2.4), desembocando en las ranuras practicadas en los largueros verticales, con lo que el único grado de libertad no impedido del ala es el giro alrededor de dicho eje; por otro lado, en su parte posterior, dos varillas atornilladas cerca del borde de salida sujetan el ala al larguero transversal trasero, pasando cada varila a través de un espárrago que fija la posición de las mismas con un tornillo de teflón (este tipo de tornillo permite fijar la posición de las varillas sin perforarlas o dañarlas). Una vez el ala es fijada dentro del túnel, para completar el montaje se debe instalar el equipo de visualización. En la Figura 2.6, donde puede apreciarse todo el montaje experimental completo, se muestra en qué consiste el equipo de visualización: la cámara de alta velocidad descrita en la Subsección 2.1.1 es instalada de forma vertical a un sistema de posicionadores que permite cambiar la posición de la cámara de forma muy precisa; además, se dispone de una lámpara de luz fría que mediante un tubo de fibra óptica permite iluminar adecuadamente la zona que se desea documentar. La presencia de la superficie libre en la zona superior del túnel hace que la visualización desde dicha zona sea imposible, por lo que debe grabarse en vertical desde abajo. 2.2 Experimentos 21 Figura 2.5 Estructura de aluminio empleada como sistema posicionador del ala dentro del túnel hidrodinámico. En su parte posterior, los largerillos incluyen dos cilindros que permiten el paso de las varillas que permiten regular la inclinación del ala. 2.2 Experimentos Tras el recorrido realizado a través del proceso de diseño y fabriación de este dispositivo generador de burbujas, es el momento de probarlo con una exhaustiva campaña experimental. En esta sección se muestran todos los detalles referentes al experimento, tales como los parámetros que se han empleado y el rango en el que estos han variado, así como el protocolo seguido para cada serie experimental. Finalmente, para concluir la sección, se especifican los resultados de las variables experimentales así como los métodos empleados para mostrar los resultados de la Sección 2.3. 2.2.1 Campaña experimental Recordando los resultados hallados en [ 17 ] y resumidos en la Sección 1.3, el diámetro y frecuencias de las burbujas producidas dependerá, dada una geometría, del caudal de gas inyectado, la velocidad del líquido exterior y del gradiente de presión local adimensionalizado en el punto de formación de la burbuja. Son estos, por lo tanto, las tres variables que configuran el espacio paramétrico que debe ser barrido en los experimentos. En la Tabla 2.3 se muestran los valores empleados de estos tres grados de libertad para los diferentes experimentos; nótese que, por las características del montaje anteriormente descrito, no se puede controlar directamente el valor del caudal, ya que este depende de la diferencia de presión entre el interior y el exterior del ala, siendo esta última variable con la velocidad y el ángulo de ataque. Por lo tanto, se toma como variable libre el valor de la presión interior del ala, debiendo calcularse el caudal de gas inyectado a posteriori. Conviene destacar los siguientes aspectos de los valores de la tabla: • Los valores de la velocidad inicidente son determinados, como ya se especificó, de forma indirecta, ya que el regulador de potencia del túnel no permite controlar directamente la velocidad. Los detalles del sistema de medición de velocidad se especifcan más adelante en la Subsección 2.2.2. • Los valores de la presión interior del depósito fueron fijados durante una de las primeras pruebas experimentales, comprobando cuándo se producían variaciones significativas en los diámetros de las burbujas al variar la presión interior del ala y manteniendo la velocidad y el ángulo de ataque constantes. • Los valores del ángulo de ataque que figuran en la Tabla 2.3 han sido calculados a posteriori según los métodos especificados más adelante. Los valores orientativos del ángulo de ataque que se pretendían conseguir varían cada 2◦ entre 4◦−12◦ . Cabe decir que no se consideraron valores α>12◦ con el fin de estar suficientemente lejos de los valores que propiciarían el desprendimiento de la capa límite. 22 Capítulo 2. Ala acuática bidimensional (a) Montaje del ala en el túnel. (b) Equipo de visualización. Figura 2.6 Montaje experimental completo. En la Figura 2.6a puede apreciarse el sistema de posicionamiento del ala y variación del ángulo de ataque mientras que en la Figura 2.6b se muestra el equipo de visualización. Tabla 2.3 Rango paramétrico empleado durante la realización de los experimentos.. Parámetros U∞[m/s]Pint [mbar]α[◦] Valores [0.3, 0.5, 0.6, 0.72] [120, 207, 288, 400, 565, 789, 890] [4.66, 6.11, 8, 10, 12] Por otro lado, para la grabación de imágenes a alta velocidad, cabe destacar que se ha empleado una resolución grabando a una velocidad de adquisición de 5000-10000 imágenes por segundo (5000-10000 fps). El siguiente paso es describir cuál es el protocolo que permite la obtención de datos de una serie experimental (esto es, una combinación de la terna [Pint ,U∞,α] ), de modo que sirva de ejemplo del proceso completo de obtención de los datos que serán procesados en la siguiente sección. El protocolo de experimentación seguido busca efectuar los realizar los ensayos de la forma más rigurosa y eficiente posible, optimizando los tiempos requeridos para cada tarea. Una vez se tiene el ala fijada en el túnel hidrodinámico, el circuito de presión conectado al ala con una presión arbitraria (pero lo suficientemente alta como para evitar que entre agua a través de los orificios) y el sistema de visualización correctamente configurado y enfocado, se procede del siguiente modo: 1. Se fija el ángulo de ataque, α , deseado para dicho experimento. Para ello se deben mover las varillas hasta encontrar (ayudándose de un transportador) aproximadamente el ángulo buscado y posteriormente 2.2 Experimentos 23 fijarlas. 2. Una vez encontrado el ángulo de ataque, se toma una fotografía de la forma más ortogonal posible de la sección transversal del ala. Dicha imagen será utilizada para calcular el ángulo de ataque real con el que incide la corriente. 3. En este punto, si fuera necesario, se vuelve al paso 1 hasta encontrar un valor aproximado al esperado para α. 4. Fijado el valor de α, se selecciona la velocidad a la que se realizará el ensayo. 5. Finalmente, se selecciona una presión interior de las que figuran en la Tabla 2.3 y se deja transcurrir un tiempo suficiente hasta alcanzar un régimen estacionario6. 6. Se configura la cámara para que inicie la grabación y se toman 3 vídeos de 3 instantes diferentes con la intención de poder realizar una media de los valores obtenidos. 7. Mientras se produce la grabación de los videos, se regula la presión para poder continuar con el siguiente valor. 8. Finalizadas todas las pruebas variando la presión interior del ala y manteniendo la velocidad y el ángulo de ataque constantes, se procede a cambiar el valor de la velocidad y repetir los pasos 5-7, y así hasta completar todas las velocidades. 9. Finalmente, si se desea, se modifica el valor del ángulo de ataque y se repiten los pasos anteriores. 2.2.2 Datos obtenidos y métodos de postproceso Llegados a este punto y concluidos todos los experimentos, se dispone del material necesario para que, tras el postproceso adecuado, se pueda realizar el análisis de los resultados. Es objeto de esta sección detallar los métodos que se han utilizado para llegar a los resultados mostrados en la Sección 2.3. Videos Fruto de los ensayos realizados se dispondrá de 3 vídeos para cada terna de parámetros (α,U∞,Pint ) . Un ejemplo ilustrativo del tipo de imágen que se obtiene de la grabación con la cámara de alta velocidad puede encontrarse en la Figura 2.7 donde se muestra la producción de burbujas en un orificio situado en el borde de ataque. Debido a que se dispone de distancias calibradas en la imagen, es imediato realizar la medición de los diámetros de las burbujas empleando un software específico como el gratuito ImageJ ® o el programa comercial de la propia cámara Phantom ® . Dado que las burbujas producidas son monodispersas, en cada vídeo se realiza la obtiene el diámetro de una burbuja característica como media de dos diámetros perpendiculares entre sí, ya que la burbuja, por efecto de la superficie del ala, puede no ser perfectamente esférica a la salida del orificio. El diámetro finalmente asignado a la terna de parámetros considerada será la media aritmética de los 3 vídeos disponibles. Cabe destacar que, dado que el proceso de medición es realizado de forma manual, se ha considerado un error de 1 pix en la medida del diámetro, lo que con la resolución y escala tomada corresponde a errores de ±25 µm . Por otro lado, para la medida de la frecuencia de burbujeo, se fija una determinada posición en la imágen y se anota el número de frames que son necesarios para que pasen por dicha posición un total de 20 burbujas. La frecuencia obtenida es por lo tanto fb=20 ×frame rate n◦de frames (2.5) Velocidad incidente Como se comentó en secciones anteriores, la velocidad no puede ser impuesta con un valor prefijado, sino que debe determinarse su valor mediante análisis de imagen. Para ello, se emplean pequeñas partículas trazaderas aproximadamente esféricas con una densidad similar a la del agua y se graban recorrer una distancia predeterminada en el túnel hidrodinámico cuando este está libre de cualquier sólido. Este proceso se repite varias veces para cada nivel de potencia del túnel seleccionads, de forma que la medida de la velocidad sea lo más robusta posible. 6 Durante la operación experimental este tiempo podía consistir en 20-30s, que es el tiempo que se tardaba en configurar y preparar la cámara para la grabación. Este tiempo es mucho mayor que el tiempo de difusión viscosa D2/νg y que el tiempo característico que emplea el gas en recorrer el tubo capilar (∼O(Lt/U∞), comprobándose esto último en la Sección 2.3 24 Capítulo 2. Ala acuática bidimensional Figura 2.7 Ejemplo de imágen extraída de uno de los videos realizados con una resolución de 1280x800 a 5000 fps. La escala del vídeo es de 25mum/px . La imagen muestra un experimento con α=6.11◦ , U∞=0.49m/syPint =289mbar. Variables de interés Una vez que se han medido los diámetros y las frecuencias de formación de las burbujas de cada terna de parámetros y se conoce del valor de la velocidad incidente, se pueden calcular otras variables de interés del problema, como es el caudal total de gas inyectado, obtenido mediante la Ecuación (1.14) Qg=πd3 b 6fb Finalmente, para poder contrastar la hipótesis de la analogía de los gradientes favorables de presión creados en entorno del borde de ataque de un perfil aerodinámico con los resultados mostrados en [ 17 ] se debe disponer del coeficiente de presión y su gradiente para los distintos casos que lo requieran. En concreto, para la definición aquí considerada del coeficiente de presión, Cp(x) = p−p∞(x) 1/2ρU2 ∞ ,con x∈Σs(2.6) resulta evidente que dicho coeficiente sólo dependerá, en cada punto y para el caso de un perfil simétrico de geometría fija, del ángulo de ataque geométrico previamente fijado; nótese que el valor del gradiente adimensional del coeficiente de presión, por su parte, dependerá (además de α ) del valor empleado en la adimensionalización de la distancia, en este caso la cuerda, c . El valor de Cp ha sido determinado por el Método de Green de Elementos de Contorno descrito en el Apéndice A. Conocido el coeficiente de presión en los puntos de integración de cada panel del perfil, resulta inmediato el cálculo de su gradiente local en cada punto que se requiera mediante derivación numérica de los valores interpolados. 2.3 Resultados y discusión Finalizados todos los ensayos y procesados todos los datos obtenidos, toca desgranar y analizar el rico abanico de fenómenos que se han observado. El objetivo de esta sección será estudiar y comprender las características del proceso de formación de burbujas en función de los parámetros ya mencionados (presión interior del ala -es decir, caudal de aire-, ángulo de ataque del ala y velocidad incidente) y de la posición del orificio sobre la superficie del ala. Dado que, como se comentó al final del Capítulo 1, existe la posibilidad de que ciertos orificios se encuentren en zonas de gradiente favorable de presión para ciertos α y en zonas de gradiente adverso para otras, se va a analizar primero el efecto de la posición del orificio de inyección en el proceso de formación de burbujas. De este modo el análisis puede centrarse en estudiar el efecto que la variación 2.3 Resultados y discusión 25 Figura 2.8 Efecto de la posición del orificio en el proceso de generación de burbujas. En la figura se muestran dos orificios situados en x/c=0.0017 a la izquierda y x/c=0 a la derecha para un ángulo de ataque, α≃12◦ , una presión interior, Pint =200mbar y 3 velocidades diferentes. Puede apreciarse que, mientras que en el orificio de la derecha la formación de burbujas permanece constante y monodispersa con independencia del valor de la velocidad, en el orificio de la izquierda se forma un ligamento gaseoso cuya extensión y tamaño aumenta a medida que lo hace la velocidad. de los tres parámetros tienen sobre los diámetros y frecuencias de las burbujas producidas centrándonos exclusivamente en un único orificio. 2.3.1 Efecto de la posición del orificio de inyección El proceso de generación de burbujas depende de la posición del orificio de inyección. En efecto, se muestra en la Figura 2.8 una imagen extraída de un ensayo realizado a ángulo de ataque α≃12◦ para tres velocidades diferentes ( U∞= [0.33,0.46,0.72]m/s ) manteniendo la presión del interior del ala constante. Como se aprecia en la figura, el orificio de la derecha, cuya posición longitudinal adimensional es x/c=0 , se mantiene en un régimen de burbujeo monodisperso independientemente del valor de la velocidad aplicado, mientras que en el orificio de la izquierda, situado en x/c=0.0017 , el comportamiento es bien distinto: para velocidades bajas, el proceso de formación de burbujas no termina de producirse cerca del orificio como sí ocurre con el orificio anterior, sino que puede apreciarse la formación de un chorro de aire eyectado desde el orificio. Este chorro crece en extensión y en diámetro aguas abajo del orificio a medida que se aumenta la velocidad. Para comprender como la formación de burbujas en orificios situados relativamente cerca el uno del otro puede presentar semejantes diferencias, se muestra en la Figura 2.9 las posiciones de los puntos de remanso y picos de succión para diferentes ángulos de ataque junto con las posiciones de los orificios de inyección. Como puede apreciarse en el diagrama, donde se ha señalado específicamente las posiciones de los orificios mostrados en la Figura 2.8, para un ángulo de ataque de α=12◦ , el orificio en la posición x/c=0.0017 se encuentra muy cerca del pico de succión calculado numéricamente con el Método de Elementos de Contorno. Debido a que la resolución del problema de un perfil bidemensional por este método constituye una buena aproximación, pero al fin y al cabo, una aproximación de la solución real del campo de velocidades y presiones alrededor del ala, parece probable que este orificio se encuentre realmente aguas abajo del pico de succión para α=12◦ y por lo tanto en una región de gradiente desfavorable de presión, por lo que la deceleración que el líquido experimenta en esta zona provoca la recirculación del gas aguas abajo del punto de inyección, formando la característica forma mostrada en la Figura 2.8. Por lo tanto, queda demostrado que la posición que ocupen los puntos de inyección de gas en el ala resulta crítica si se quiere asegurar una generación de burbujas exitosa, ya que una pequeña desviación puede hacer que un orificio, para determinados ángulos de ataque, pase de un régimen de burbujeo monodisperso a un régimen en el que forme grandes masas de aire. Cabe entonces preguntarse en este punto cómo varía el coeficiente de presión y más concretamente su gradiente para orificios a lo largo de la coordenada curvilínea con origen en el borde de ataque del perfil y 32 Capítulo 3. Ala acuática rotatoria Tabla 3.1 Medidas del tanque empleado en los experimentos del ala acuática rotatoria. Dimensión Longitud [m] Altura [m] Espesor [mm] Valor 1 0.6 30 • Además debe existir una distancia adecuada a las paredes del tanque con el fin de evitar impactos y velocidades inducidas por el mismo. • Por otro lado, las palas deben estar formadas en su sección transversal por perfiles aerodinámicos donde existan fuertes gradientes favoralbes de presión entre el punto de remanso y el pico de succión de cada sección transversal para cada ángulo de ataque. • Con el fin de explorar diferentes resultados en función del gradiente de presión, el ángulo de ataque debe poder ser regulable. • Finalmente, se debe conseguir suministrar aire desde el interior de las alas hacia el líquido, lo que deberá hacerse a través del eje de rotación, convirtiéndo al igual que en el caso bidimensional cada pala en un depósito. Para cumpliar todas estas premisas, se va a describir el diseño de cada uno de los componentes por separado. 3.1.2 Banco de ensayos El banco de ensayos para un dispositivo de las características mencionadas arriba debe consistir fundamentalmente en un tanque lleno de líquido, en este caso nuevamente agua tanto por su sencillez como por la futura aplicación del propio dispositivo 1 . Por otro lado, el tanque debe realizarse en un material transparente que permita la visualización de las burbujas producidas con el fin de poder medir los diámetros y las frecuencias de producción como se realizó en el Capítulo 2; se ha optado en este caso por el metacrilato como material de fabriación, debido a su menor coste con respecto al cristal. Aunque un tanque de geometría cilíndrica sería una solución idónea teniendo en cuenta la axilsimetría del problema, la distorsión óptica provocada por la combinación de la geometría curva y el espesor del metacrilato nos llevan a elegir un tanque de sección cuadrada. En cuanto a la altura, se ha considerado suficiente dipsponer de una altura total del tanque de 60 cm, con lo que la altura total de agua estará en torno a los 0.5 m. Con estos requerimientos, y basándose en la experiencia del fabricante, el espesor de pared recomendado no debería ser inferior al mostrado en la Tabla 3.1, con el fin de aguantar los casi 500 kg del peso del agua. Por otro lado, dado que el equipo de visualización empleado será el mismo y una vez más la interfase aire-agua imposibilitaría la visión desde arriba, la grabacion de imágenes a alta velocidad ha de efectuarse desde la zona inferior del tanque, por lo que este debe estar elevado. Para poder realizar los experimentos se requiere por tanto idear una estructura que pueda ser utilizada como banco de ensayos, debiendo la misma dar soporte para la realización de las siguientes tareas: • Soporte del tanque a una altura suficiente como para poder incluir el montaje completo del equipo de visualización bajo el mismo. Esta altura debe ser la mínima posible con el fin de evitar aumentar la longitud de los pilares verticales y con ello la inestabilidad por pandeo. • Soporte para el equipo rotativo en la zona superior del tanque, el cual debe soportar el peso del motor y sistema de aireación junto con los esfuerzos radiales producidos por el giro del motor. • Estabilidad ante esfuerzos en la dirección perpendicular del eje. En efecto, el movimiento creado por las palas aerodinámicas en el seno del líquido harán que el tanque experimente esfuerzos en su pared que serán, en última instancia, transmitidos a la estructura que constituye el banco de ensayos. De este modo, teniendo en cuenta todo lo anterior, se ha diseñado y montado una estructura con perfiles de aluminio, que aporta la rigidez y robustez necesarias para satisfacer los requerimientos arriba descritos. En la Figura 3.1 se muestra una perspectiva de la estructura y el tanque aquí descritos, donde puede observarse que el tanque se sitúa aproximadamente a 1 m de distancia del suelo, dejando espacio suficiente para el montaje del equipo de visualización. Además los dos últimos estantes proveen soporte para el motor eléctrico y el depósito de suministro de presión (este útlimo descrito más adelante), mientras que los perfiles en diagonal aportan la estabilidad suficiente para que, en operación, la estructura no se alabee. 1 Conviene mencionar que el empleo de agua corriente del grifo, a pesar de no ser agua pura, no tiene nada que ver con el tipo de aguas encontrado en la industria de tratamiento de aguas, debido a la abundancia de partículas y altos contenidos de ácidos como nitratos. 3.1 Diseño y fabriación 33 Figura 3.1 Perspectiva del banco de ensayos empleado en la campaña experimental. En las figura se muestra la estructura fabricada a base de perfiles de aluminio junto con el tanque de metacrilato en su interior. 3.1.3 Diseño de las alas El diseño de las alas empleadas para este dispositivo será muy similar al empleado en el dispositivo bidimensional, por lo que muchos de los criterios de diseño seguidos en la Subsección 2.1.2 pueden ser extrapolados a este caso. Así, por ejemplo, el perfil aerodinámico empleado será el mismo que en el caso bidimensional, es decir el perfil simétrico NACA 0012. Sin embargo, sí que existen una serie de diferencias que merece la pena pararse a detallar. La primera de ellas es que ahora no se dispone de restricciones especiales ni para la cuerda ni para la envergadura del ala, ya que cada una debe ir unida al eje a través de su brazo correspondiente, existiendo por lo tanto un único punto de sujecion situado cerca de la zona de mayor espesor de la sección transversal del ala. En este caso, dado el carácter tridimensional del problema y con la intención de poder aplicar la teoría aerodinámica general para alas esbeltas (descrita en detalle en el Apéndice B ), sería conveniente que su alargamiento, definido como Λ=b2/S con b la envergadura y S la superficie en planta, sea tal que Λ>1 . Para el caso de un ala rectangular, la expresión del alargamiento es simplemente el cociente entre la envergadura y la cuerda, Λ=b/c , por lo que para el diseño de las palas considerado, una de las restricciones será que b>c. Por otro lado, la sección transversal del ala será muy similar a la del caso bidimensional del Capítulo 2, ya que el ala actuará como depósito respecto a los orificios de inyección. No obstante, en este caso la inyección de aire se realiza desde un lateral del ala en lugar de hacerla desde la zona del extradós, aprovechando el punto de conexión del ala al eje de rotación a través del brazo. En cuanto a los orificios de inyección, estos consistirán al igual que en el caso del ala bidimensional en idénticos tubos capilares de acero con ratios D/Lt1 , con el fin de conseguir que la perdida de carga evite variaciones bruscas de caudal por el cambio en la velocidad o el ángulo de ataque incidentes. Una importante diferencia con respecto al caso del perfil bidimensional, es que en este caso no se explorarán diferentes zonas de inyección en torno al borde de ataque, sino que empleando los resultados mostrados en la Subsección 2.3.1, se ha decidido colocar los puntos de inyección justo en el borde de atque de cada pala. Además, se ha colocado un número mayor de puntos de inyección con la intención de analizar el comportamiento en distintos puntos de la envergadura y de emular aún más el posible de diseño de un prototipo de dispositivo de aireación y agitación de una planta depuradora. 34 Capítulo 3. Ala acuática rotatoria Tabla 3.2 Dimensiones características relevantes del diseño del ala y los orificios de inyección. El valor de g en la última columna representa la distancia de separación entre los diferentes orificios a lo largo de la envergadura. Elemento b[mm]c[mm]D[µm]Lt[mm]gap[mm] Valor 200 125 160 30 5 (a) Perspectiva en CAD del ala, con un corte longitudinal. (b) Perspectiva final. Figura 3.2 Diferentes perspectivas del diseño y fabricación del ala para el prototipo de dispositivo rotatorio. Finalmente, el diseño del ala queda completo especificando el sistema de variación de ángulo de ataque de las palas. Dado que existe un único punto de sujeción y por este pasarán no sólo los esfuerzos causados por la fuerza centrífuga sino también los ejercidos por la presión interior dentro de los tubos que confeccionan los brazos del dispositivo, parece conveniente que la unión entre las alas y los citados brazos metálicos no se haga directamente atornillando sobre la pala, puesto que el material de esta continúa siendo un polímero de ABS. En su lugar, se ha confeccionado una solución consistente en la fabricación de un perfil realizado en aluminio idéntico al del ala, pero con un espesor mucho menor que la envergadura de esta. De este modo, el perfil sirve de nexo de unión entre el ala y eje de rotación a través del brazo. Esta pieza de aluminio, tal y como se observa en la Figura 3.2 posee un total de 4 orificios avellanados, 3 de los cuales permiten una unión robusta del ala y el otro se atornilla mediante un tornillo hueco directamente al brazo, donde se ha realizado una rosca interna para tal fin. De este modo, simplemente girando el perfil de aluminio y apretando el mencionado tornillo, se tiene un sistema capaz de implementar el ángulo de ataque de la pala deseado y de, a su vez, permitir el paso de aire a través de sí. Téngase en cuenta, además, que el ángulo de ataque no puede ser prefijado, sino que debe medirse de forma adecuada mediante análisis posterior de imágen, similarmente al caso bidimensional. 3.1.4 Equipo rotativo y suministro de presión Finalmente, se llega al último punto a completar en el diseño del prototipo dispositivo agitador y generador de microburbujas: el sistema de rotación y aireación. El sistema de rotación es una tarea sencilla de implementar, pues basta un motor eléctrico y un eje lo sucicientemente largo (aunque no demasiado si se quieren evitar excentricidades de la carga); la única condición "especial" que debe satisfacer dicho eje es que consista en un tubo hueco por dentro, de forma que permita el paso de aire desde el circuito de presión hasta las alas. El motor, por su parte, posee un regulador de frecuencia que permite controlar la potencia suministrada e indirectamente velocidad de rotación; la velocidad de rotación tendrá que medirse por lo tanto a posteriori. El diseño de un adecuado sistema de suministro de presión, por otro lado, es una tarea un poco más compleja. El sistema que se ha propuesto para permitir al ala rotar al tiempo que se transporta el aire desde la red de presión hasta el interior del ala consiste en lo siguiente: • Un recipiente de aluminio fabricado específcamente para el dispositivo servirá de depósito estanco, situado fuera del tanque. • El depósito se encontrará alineado con el eje del motor, de forma que éste último pueda pasar a través de él. • En la tapa superior e inferior del depósito se disponen 2 rodamientos (uno a cada lado) junto con retenes para permitir al eje rotar manteniendo el depósito estanco en todo momento. • El eje posee en la zona que se encuentra en el interior del depósito unas ranuras que permiten la entrada de aire desde el depósito presurizado hacia su interior. Este aspecto implica además que el 3.2 Experimentos 35 (a) Corte transversal del interior del depósito. (b) Perspectiva final del depósito. Figura 3.3 Ilustración del diseño y fabricación del depósito de rotación estanca implementado en el dispositivo generador de microburbujas propuesto.. eje permanece siempre fijo respecto del depósito, esto es, tiene impedidos todos los desplazamientos incluido el paralelo a sí mismo. • El depósito se conecta directamente a un manorreductor conectado a su vez a la red de presión, por lo que la presión interior del depósito y por lo tanto del ala es controlable desde este manorreductor de forma similar a como se hizo en el dispositivo del Capítulo 2. En la Figura 3.3 pueden observarse algunos de los detalles del diseño del depósito y de su fabriación así como el aspecto de la configuración final con el motor eléctrico. 3.1.5 Ensamblaje final Finalmente, descritos todos los componentes que conforman el banco de ensayos y el prototipo de agitación y aireación con microburbujas que se propone, se muestra en la Figura 3.4 el montaje experimental completo, donde se ha incluido además el equipo de visualización situado bajo el tanque, listo para comenzar la campaña experimental. 3.2 Experimentos Una vez completado el diseño de todo el prototipo se está en disposición de comenzar la campaña experimental que proporcione los datos necesarios para evaluar la validez de las leyes de escala obtenidas en el capítulo anterior. De este modo, esta sección se estructura de forma análoga a su homóloga en el Capítulo 2, mostrando en primer lugar el espacio paramétrico que se pretende explorar durante la campaña experimental y posteriormente describiendo los métodos de análisis y postproceso utilizados. Cabe destacar que estos métodos supondrán una ligera modificación/ampliación de los ya expuestos en el capítulo anterior. 3.2.1 Campaña experimental Las variables de las que dependen el diámetro y las frecuencias de producción de las burbujas para el problema tridimensional son las mismas que para el problema bidimensional, por lo que el el tipo de parámetros y el rango explorado será muy similar al del problema 2D. Las variables libres que en este problema se pueden controlar son las siguientes: • Velocidad de rotación de las palas. El valor que se controla es un parámetro arbitrario de frecuencia, por lo que la velocidad de rotación debe extraerse de las imágenes grabadas. • Presión en el interior del depósito. Una vez más, no es posible el control directo del caudal de inyección de aire sino que sólo se puede actuar sobre la presión interior del depósito; no obstante, la pérdida de carga generada por los tubos capilares implicará, como en el caso bidimensional, que las variaciones 36 Capítulo 3. Ala acuática rotatoria Figura 3.4 Perspectiva general del montaje experimental completo. Tabla 3.3 Espacio paramétrico a explorar en la campaña experimental. La última columna representa las posiciones de los orificios que son analizados en cada vídeo, donde rt representa la distancia de cada orificio al eje de rotación, por lo que la velocidad incidente para una velocidad de rotación dada, Ωi, y para un orificio concreto, rtjes Vi j =Ωirtj. Parámetros Ω[rpm] Pint [mbar] αg[◦]rt[mm] Valores [65,80,105,125] [125,150,175,200,250,300] [8,12] [389,384,379], [424,419,414] ,[469,465,460] de presión provocadas por el aumento de velocidad del líquido y/o por el cambio de ángulo de ataque originen variaciones pequeñas del caudal de aire. • Ángulo de ataque geométrico, αg . Este es el ángulo de ataque que poseen las palas con respecto a un plano de normal al eje de rotación. Aunque en este caso todos los orificios han sido situados en el borde de ataque del ala, los resultados obtenidos para el diámetro y las frecuencias variarán de un orificio a otro. En efecto, si bien la presión en el interior del ala puede considerarse uniforme y el ángulo de ataque geométrico es constante para todas las secciones de la pala, existen dos magnitudes que varían con la envergadura: la velocidad local de cada sección (producto de la velocidad angular de giro, Ω , y de la distancia de cada orificio al eje de rotación, r , y el ángulo de ataque efectivo, αe f f , definido en la teoría general de alas como αe f f =αg−αind (y) , siendo αind (y) el ángulo de ataque inducido por los torbellinos de punta de pala que varía a lo largo de la coordenada que recorre la envergadura, y , en el sistema de ejes viento tal y como se detalla en el Apéndice B. Por lo tanto, parece necesario, para cada terna de parámetros ( Ω,α,Pint ) se analicen los diámetros y las frecuencias de distintos orificios a lo largo de la envergadura con el fin de poder capturar toda la fenomenología implicada en el problema. El espacio paramétrico explorado durante la campaña experimental queda resumido en la Tabla 3.3 El protocolo experimental seguido en este caso para obtener una serie experimental es el siguiente: 3.2 Experimentos 37 Figura 3.5 Muestra de imágen extraídas de un vídeo realizado con αg=8◦ , Pint =150mbar y Ω=80rpm . La resolución original del vídeo para este caso es de 688x456 con una escala de 21 µm/px grabado a 20.000 fps. 1. Selección del ángulo de ataque geométrico a ensayar. El ángulo de ataque geométrico de cada pala debe ser opuesto con respecto al plano de normal según el eje de rotación, ya que se pretende evitar momentos generados en el eje por efecto de la sustentación. 2. Verificación mediante análisis de imágen del ángulo de ataque geométrico. 3. Colocación de las palas e inmersión en el tanque de agua, conectando el sistema de presión para evitar la entrada de agua en el interior de las alas. 4. Posicionamiento de la cámara enfocando alguno de los 3 rangos de orificios especificados en la Tabla 3.3 5. Selección de la menor velocidad de rotación a analizar. 6. Selección de la presión en el interior del depósito. 7. Accionamiento del mecanismo de rotación y espera durante el tiempo requerido para que el sistema se considere estacionario2. 8. Captura vídeo donde, al menos, haya transcurrido más de una vuelta completa del ala. Esto se debe a que el ala sobre la que se visualizan los resultados es aquella que posee un ángulo de ataque geométrico negativo con respecto al plano de normal según la dirección del eje de rotación y sentido hacia el motor, por lo que se requiere al menos una vuelta completa para poder efectuar la medida de la velocidad de rotación. 9. Una vez terminado, cambio en la presión interna del depósito con la misma velocidad. 10. Una vez recorrido todo el rango de presiones, cambio en la velocidad de rotación. 11. Finalizadas todas las velocidades, volver al paso 3 y repetir todos los pasos hasta completar todos los orificios de interés. 12. Finalmente, cambiar el ángulo de ataque geométrico y repetir el proceso. 3.2.2 Datos obtenidos y postproceso Una vez se haya finalizado la campaña experimental se dispondrán de 3 vídeos para cada orificio y terna de parámetros ( Ω,αg,Pint ). Un ejemplo ilustrativo del aspecto de estos vídeos puede observarse en la Figura 3.5 donde puede comprobarse que, al igual que en el caso bidimensional, se puede emplear software de tratamiento de imágenes para medir de forma manual el diámetro de las burbujas y las frecuencias de producción, exactamente de la misma forma que se realizó en la Subsección 2.2.2 y asignando un error de 1 px en la medida del diámetro, que se corresponde en este caso teniendo en cuenta la resolución de la Figura 3.5 con ±21 µm/px . La medida de la velocidad angular de giro, por otro lado, puede obtenerse del número de frames transcurridos entre un paso de la pala y el siguiente. La medida del caudal en el problema tridimensional se efectúa de forma indirecta del mismo modo que en el caso bidimensional, esto es, aplicando la Ecuación (1.14) que relaciona diámetros de burbuja y frecuencia de formación de burbujeo con el caudal. Hasta este punto el postprocesado de resultados resulta completamente análogo al descrito para el caso bidimensional, con la salvedad de que el número de orificios donde hay que realizar las medidas es mayor. La única variable que resta por calcular ahora es el gradiente de presión local en el orificio, es decir Ps(¯s=0) . Sin embargo, no es posible emplear el código de elementos de contorno del Apéndice A directamente utilizando 2 Debe apreciarse en este caso el caracter claramente no estacionario del problema, por lo que se asume que, transcurrido un tiempo en el que todos los parámetros mantienen un valor constante, se alcanza un régimen permanente. 38 Capítulo 3. Ala acuática rotatoria el ángulo de ataque geométrico, ya que para el caso de un ala de envergadura finita el ángulo de ataque efectivo en cada sección del perfil se ve reducido debido al wash-down creado por los torbellinos de punta de pala que inducen velocidades verticales a lo largo de toda la envergadura del ala. Por ello, el ángulo de ataque efectivo con el que la corriente incide en cada sección se puede modelar como αe f f (y) = αg−αind (y) para el caso de un ala sin torsión. El ángulo de ataque inducido, αind (y) , es pues una medida del ángulo de ataque que se deflecta la corriente incidente debido a los torbellinos de borde marginal. En efecto, en esta región, la presión en el extrados e intradós se iguala, provocando el movimiento del fluido desde el entorno del borde marginal del intradós hacia el extradós, generando así un torbellino en la punta. Debido al caracter finito de la envergadura, este torbellino induce velocidades transversales a la corriente incidente dirigidas hacia la raíz sobre la superficie del ala, al tiempo que deflecta la corriente incidente un ángulo αind =wi/U∞ , con wi la velocidad vertical que crea el torbellino de borde marginal y que puede escribirse de la forma wi(y) = 1 4πZy=b/2 y=−b/2 dΓ/dy0 y−y0 dy0(3.1) siendo Γ(y) = Hv·dl el valor de la circulación alrededor del perfil situado a una distancia y de la raíz del ala. La deducción de la Ecuación (3.1) junto con las definiciones, hipótesis realizadas y otra serie de razonamientos de obligada comprensión sobre la aerodinámica de alas de envergadura finita puede encontrarse en el Apéndice B. El ángulo de ataque inducido para cada posición y de la envergadura puede obtenerse como αind (y) = 1 4πU∞Zy=b/2 y=−b/2 dΓ/dy0 y−y0 dy0(3.2) En el Apéndice B se incluyen los detalles de un método numérico conocido como Vortex-Lattice que permite resolver problema sustentador de un ala de envergadura finita. La resolución de este problema PROPORCIONAel valor de la circulación total Γ , de la distribución de circulación, Γ(y) y del ángulo de ataque inducido, αind (y). De este modo, para el caso de un ala aislada con ángulo de ataque geométrico αg sobre la que incide una corriente uniforme U∞ , se podría emplear el código de elementos de contorno empleado en el Capítulo 2 suministrando al programa el valor del ángulo de ataque efectivo para el perfil de la sección que se desee analizar. Sin embargo, en este problema, existen dos diferencias con respecto al caso de un ala sobre la que incide una corriente uniforme a ángulo de ataque geométrico αg : la primera es que en nuestro problema la velocidad varía de forma lineal con la distancia al eje y la segunda es que la velocidad que incide sobre una pala se encuentra perturbada por la estela que deja tras de sí la otra pala. Para solventar el primer inconveniente se propone un modelo en el que se considera la corriente incide sobre el ala con una velocidad uniforme de valor U∞≃Ω(rmax +rmin)/2 , es decir, una velocidad media entre la máxima y mínima a lo largo de la envergadura. Por otro lado, el modelado de la estela tras el ala no es un problema trivial aún hoy en día. Dado que un modelo que reproduzca de forma rigurosa el flujo en la estela queda fuera del alcance de este trabajo, se va a proponer uno más sencillo que nos permita obtener el ángulo de ataque efectivo del ala. Considérese la Figura 3.6, donde se representan de forma esquemática las dos palas y las estelas generadas tras de sí. Como se ha comentado, los torbellinos de borde marginal inducen velocidades transversales sobre la superficie del ala y verticales sobre la corriente incidente que modifican el ángulo de ataque efectivo de cada sección de la pala. Estos torbellinos se convectan aguas abajo en la estela, de forma que los torbellinos del ala que va por delante provocarán una reducción del ángulo de ataque efectivo del ala de atrás del orden de ∼αind (y) , dependiendo esta disminución de la sección, y , considerada. Así pues, dado que el ángulo de ataque efectivo para un ala que enfrenta un flujo uniforme es αe f f (y) = αg−αind (y), se va a modelar el efecto de la estela considerando que el ángulo de ataque efectivo total del ala es αe f f (y) = αg−2αind (y)(3.3) 3.3 Resultados y discusión En esta última sección se van a analizar los resultados de la campaña experimental descrita en la sección anterior utilizando para ello los métodos anteriormente descritos. La estructura de la sección será similar a la que se siguió en el Capítulo 2: en primer lugar, se presentan de forma cualitativa los resultados obtenidos a través de los vídeos tomados para distintos valores representativos de los parámetros del problema de forma 3.3 Resultados y discusión 39 Figura 3.6 Representación esquemática del problema de determinación del ángulo de ataque efectivo de las palas. En la figura se muestra como los torbellinos de borde marginal se convectan aguas abajo en la estela. La estela de un ala modifica el flujo incidente sobre la otra, reduciendo aún más el ángulo de ataque efectivo del ala posterior. (a) Borde marginal de la pala. (b) Zona central de la pala. (c) Raíz de la pala. Figura 3.7 Resultados obtenidos en 3 zonas diferentes a lo largo de la envergadura para unas palas rotando a un ángulo de ataque α=12◦yΩ=80rpm.. que, tras la cuatificación de los resultados, se pueda comprobar si un escalado similar al realizado al caso bidimensional es también válido para este problema. 3.3.1 Fenomenología Esta sección se centrará en mostrar cualitativamente los fenómenos representativos hallados durante la campaña experimental. En primer lugar, a diferencia del caso bidimensional, cabe destacar que los orificios están todos situados en la punta del borde de ataque de los perfiles que conforman las alas. Sin embargo, dado el caracter finito de la envergadura de las alas y dado que estas se encuentran rotando respecto a un eje, sí deben existir diferencias entre los distintos puntos de inyección a lo largo de la envergadura. La Figura 3.7 muestra de forma cualitativa las diferencias existentes para distintos orificios distribuidos a lo largo de la envergadura de la pala para un ángulo de ataque geométrico α=12◦ , una velocidad angular Ω=80rpm y una presión interior Pint =150mbar . En principio, cabría esperar que dado que la velocidad en la punta de la pala es mayor que en la raíz debido a la rotación, las frecuencias de burbujeo para un caudal constante deberían aumentar y consiguientemente los diámetros de las burbujas formadas disminuir. No obstante, como se observa en las diferentes imágenes, no se aprecian a simple vista diferencias demasiado significativas. Analizada la influencia de la posición de la envergadura la Figura 3.8 muestra como cambian los diámetros de las burbujas al variar el caudal de gas inyectado. Como cabía esperarse, al aumentar el caudal de aire inyectado el diámetro de las burbujas aumenta, algo análogo a lo visualizado en el caso bidimensional. Finalmente, dada la influencia que ejerce el gradiente de presión local sobre los diámetros y frecuencias de las burbujas para el caso bidimensional, cabe preguntarse si en este caso se obtendrán resultados similares a los hallados en la Sección 2.2. La Figura 3.9 muestra los orificios de la zona central del ala para los dos ángulos de ataque considerados durante la campaña experimental (8 y 12 grados). Puede apreciarse que la influencia del ángulo de ataque geométrico es distinta al caso de un perfil bidimensional, donde las reducciones en los diámetros por efecto de aumento del gradiente de presión al variar el ángulo de ataque podían llegar a ser del orden de 150 micras. 40 Capítulo 3. Ala acuática rotatoria (a) Pint =125 mbar. (b) Pint =250 mbar. Figura 3.8 Influencia del caudal de gas inyectado en los diámetros de la burbuja. Las figuras muestran los las burbujas producidads en orificios de la zona central del ala para un α=8◦yΩ=65rpm. (a) α=8◦.(b) α=12◦. Figura 3.9 Influencia del ángulo de ataque geométrico en los diámetros de las burbujas producidas. Ambas figuras muestran el caso de orificios situados en la zona central del ala cuando la presión del depósito es Pint =150mbar y la velocidad de rotación Ω=80rpm. . 3.3.2 Resultados cuantitativos Una vez visualizada la fenomenología encontrada durante la realización de los experimentos es momento de evaluar cuantitativamente los datos resultados obtenidos. Antes de comenzar a visualizar los valores de los diámetros de burbujas obtenidos para este prototipo, conviene cuantificar el valor del ángulo de ataque efectivo para cada sección del perfil. La Figura 3.10 muestra el valor del ángulo de ataque efectivo, obtenido tal y como se especificó en la Subsección 3.2.2, esto es, restando al ángulo de ataque geométrico el doble del ángulo de ataque inducido calculado numéricamente en cada sección del ala. Este sencilo modelo pretende tener en cuenta en el cálculo del gradiente local de presión el efecto de los tobellinos de borde marginal presentes en alas de envergadura finita y el efecto de la estela que el otro ala ejerce sobre el flujo incidente (situado con ángulo de ataque opuesto respecto a su radialmente simétrica). Como puede comprobarse, el ángulo de ataque efectivo es mucho menor que el ángulo de ataque geométrico que se pretendía imponer y, lo que es más importante, debido al incremento en el ángulo de ataque inducido cuando se aumental el ángulo de ataque geométrico, las diferencias en el ángulo de ataque efectivo en cada sección para los valores de los ángulos geométricos considerados son prácticamente despreciables. Este hecho explicaría porqué no se apreciaban variaciones significativas en los diámetros de las burbujas cuando se observaban a dos ángulos de ataque muy diferentes (8 y 12 grados). Ahora se está en condiciones de analizar cuantitativamente las variaciones en los diámetros de las burbujas que existen a lo largo de la envergadura. La Figura 3.11 muestra los diámetros de burbuja medidos para los diferentes puntos de inyección considerados durante los experimentos y para los dos ángulos de ataque analizados. En la Subsección 3.3.1 se comentó que, en principio, cabría esperarse encontrar una reducción en el diámetro de la burbuja a medida que nos alejamos de la raíz, debido al incremento de velocidad incidente. Sin embargo, como se observa en la Figura 3.11, en la punta se encuentran precisamente diámetros de burbuja 3.3 Resultados y discusión 41 -1 -0.5 0 0.5 1 0 1 2 3 4 Figura 3.10 Ángulo de ataque efectivo para distintos puntos a lo largo de la envergadura. En la figura se muestran los valores del ángulo de ataque efectivo cuando se le sustrae al ángulo de ataque geométrico el doble del valor del ángulo de ataque inducido calculados numéricamente utilizando el método de Vortex-Lattice descrito en el Apéndice B. mayores que en el centro de la pala, además de ser en ciertos casos más pequeños los diámetros asociados a ángulos de ataque geométrico de 8 grados que los asociados a αg=12◦ . Este aspecto encaja perfectamente con la disminución del ángulo de ataque efectivo descrita anteriormente, pues precisamente en las puntas (ambos bordes marginales en este caso) es donde la reducción por efecto del ángulo de ataque inducido es mayor debido al efecto del torbellino de punta de pala. Así pues, esta disminunción en el ángulo de ataque efectivo hace que sea en el centro de la pala donde se alcanzan menores diámetros a pesar de ser menor la velocidad que en la punta, lo que pone de manifiesto una vez más que el fenómeno de formación de burbujas está dominado por el gradiente longitudinal de presión local. Finalmente, si se adimensionalizan los resultados de la misma forma que en la Sección 2.3 y en [ 17 ] y se hace uso de las leyes de escala deducidas previamente para el caso bidimensional, esto es db∝ Qg U∞L2pPs/2!2/5 y se tiene en cuenta que el valor del gradiente adimensional de presión en el punto de inyección de gas Ps=dCp/d ¯s(¯s=0)c/L , debe ser calculado para cada orificio empleando el mismo Método de Elementos de Contorno que en el caso bidimensional pero utilizando el ángulo de ataque efectivo para cada punto de la envergadura, así como el uso de un valor promedio de la velocidad incidente, dado por Ω(rmax +rmin)/2 , se obtienen los resultados de la Figura 3.12. Cabe destacar que la expresión de la recta de mejor ajuste encontrada para este caso es db L=1· q pPs/2!2/5 (3.4) que es prácticamente idéntica a la encontrada en el caso bidimensional (Sección 2.3) y en [17]. 48 Capítulo A. Introducción a la Aerodinámica potencial característica, en el caso que nos ocupa Fr ∼O(1) , ya que las velocidades del líquido son U∞∼∞m/s . Sin embargo, la única fuerza másica aquí considerada es la gravedad, que es una fuerza conservativa, por lo que deriva de un potencial, fm=−∇U y puede incluirse dentro del término de presión de la Ecuación (A.7) , −∇(p+ρU), denominándose presión reducida y sin que ello suponga pérdida de generalidad. Otra importante simplificación es posible si se considera que en los perfiles aerodinámicos como los que se emplean en los dispositivos mostrados poseen una cuerda tal que c∼O(0.5m) , y enfrentan una corriente incidente de velocidad U∞∼O(1m/s) , con lo que se tiene que el número de Reynolds, Re , que es el número adimensional que surge de comparar el término de inercia ρv·∇v con el término de esfuerzos viscosos, n·¯ ¯ τ τ τ000 , en la Ecuación (A.7), puede estimarse como Re =ρU∞c µ∼103·1·0.5 10−3∼5·1051(A.12) por lo que se concluye que puede despreciarse, en primera aproximación, los efectos de la viscosidad a distancias del orden de la longitud de la cuerda y mayores en el dominio fluido. Por otro lado, si en la Ecuación (A.8) se compara el flujo de calor por convección, ρcvv·∇T , con el flujo de calor por conducción, k∇2T , y se define además el número de Prandtl como el cociente entre la difusión viscosa y la difusión términa, Pr =µcp/k , con cp=cv para un líquido de propiedades constantes como el aquí considerado, se tiene que ρcvv·∇T k∇2T∼ρU∞cv∆T/c k∆T/c2∼µcvρU∞c kµ∼RePr =Pe 1(A.13) donde se ha definido el número de Peclet, Pe =RePr , y se ha tenido en cuenta que Pr ∼O(7) para el caso del agua. Este resultado nos lleva a que también pueden despreciarse en primera aproximación los términos de conductividad términa frente a los de convección en la ecuación de la energía. Las estimaciones anteriores son válidas cuando se toman los gradientes a distancias del orden c del objeto, lo que permite despreciar ĺos términos de esfuerzos viscosos y de conducción térmica con errores del orden de Re−1 y Pe−1 , respectivamente. Nótese, sin embargo, que dichos términos no pueden ser despreciados en las regiones cercanas a la pared del sólido, dado que el fluido debe adaptarse para cumplir las condiciones de contorno de velocidad relativa nula e igualdad de temperaturas y de flujos de calor. Para comprobar que, efectivamente, esta región es de un tamaño muy inferior al dominio, se debe considerar que en el espesor de esta región, δ , el fluido modificará su velocidad desde la corriente exterior (del orden de U∞ ) hasta hacerse nula en la pared del sólido (siempre en los ejes ligados al cuerpo), por lo que en esta región, tanto los términos viscosos como los de inercia son comparables, ρv·∇v µ∇2v∼1→∼ ρU2 ∞/c µU∞/δ2∼1⇒δ c∼Re−1/21(A.14) por lo que, para el caso que nos ocupa, la región donde los efectos viscosos no pueden ser despreciados posee un espesor muy inferior a las dimensiones características del ala; a esta región se le conoce con el nombre de capa límite. De un modo completamente análogo, puede obtenerse que la región donde los efectos de conducción no pueden despreciarse, será una región de espesor δT c∼Re−1/2Pr−1/21(A.15) denominada capa límite térmica, también muy inferior a c. Que el espesor relativo de la capa límite sea tal que δ/c1 permite realizar una gran simplificación en la Ecuación (A.7) , de forma que junto con la Ecuación (A.9) permitan abordar con relativa sencillez el problema de la capa límite. En efecto, llamando s a la coordenada curvilínea tangente a la superficie del perfil e y a la componente ortogonal a la superficie del sólido (véase la Figura A.2), los gradientes transversales (∇T) y longitudinales (∇c) de presión han de ser, respectivamente, del orden de ∆Tp∼ρ(∆V)2,∆cp∼ρU∞∆U∞(A.16) siendo V la velocidad transversal del fluido a distancias del orden δ del objeto. Estas velocidades pueden ser estimadas empleando la ecuación de continuidad ∂u ∂s+∂v ∂y=0⇔∂v ∂y∼∂u ∂s⇒∆V∼δ c∆U∞(A.17) A.3 Ecuaciones simplificadas 49 quedando entonces que ∆Tp ∆cp∼δ c2 ∼Re−11(A.18) lo que expresa que, con errores relativos del orden de Re−1 , las variaciones de presión transversal pueden ser despreciadas en la región de la capa límite, cumpliéndose las leyes de fluidoestática que surgen de resolver la ecuación ∂p ∂y=0(A.19) Para facilitar la comprensión cualitativa del fenómeno de capa límite, se va a considerar el caso estacionario de un líquido incompresible para el que las ecuaciones de continuidad y cantidad de movimiento se reducen a ∂u ∂s+∂v ∂y=0(A.20) ρu∂u ∂s+v∂u ∂y=−∂p ∂s+µ∂2u ∂y2(A.21) donde en la segunda ecuación se ha considerado despreciable el término ∂2u/∂s2 frente a ∂2u/∂y2 . Ahora bien, dado que las variaciones transversales de presión en la capa límite son muy pequeñas comparadas con las longitudinales, en cada punto de la capa límite, la presión tendrá el mismo valor en todo su espesor que el que posee en la zona exterior a la misma, esto es en y→∞ , donde la presión puede ser estimada teniendo en cuenta que esta región los efectos viscosos son despreciables y por lo tanto es de aplicación la ecuación de Euler-Bernoulli p∞+1 2ρU2 ∞=pe(s)+ 1 2ρu2 e(s)(A.22) con pe(s) y ue(s) la presión y la velocidad tangencial al perfil en la frontera de la capa límite viscosa. Nótese ademaś que, tanto en la ecuación de Euler-Bernoulli como en la ecuación de capa límte se ha despreciado el término de fuerzas másicas, sin embargo esto no supone una pérdida de generalidad pues, tanto si se considera el caso de gases donde los efectos de inercia son mucho mayores que los de gravedad (número de Froude, Fr =U2 ∞/gc 1 )) como si se considera que las fuerzas másicas derivan de un potencial, fm=−∇U , las ecuaciones obtenidas son, en esencia, las mismas. Nótese que, gracias a la ecuación de Euler-Bernoulli, el término del gradiente longitudinal de presión puede ser descrito como −∂pe(s) ∂s=ρue(s)due ds(A.23) por lo que la ecuación de cantidad de movimiento finalmente queda ρu∂u ∂s+v∂u ∂y=ρue(s)due ds+µ∂2u ∂y2(A.24) sujeta a las condiciones de contorno y=0,u=v=0,y→∞,u→ue(s)(A.25) donde se ha considerado que, como δ/c1 , la velocidad en y→∞ debe ser muy próximma al valor de la velocidad en el exterior de la capa límite. Una simple inspección de la Ecuación (A.24) revela que el campo de velocidades en el interior de la capa límite está controlado por el valor de due/ds . Así, por ejemplo, en los casos en los que la corriente se acelere a lo largo de la coordenada s , el gradiente favorable de presión hará que las partículas en el interior de la capa límite se aceleren a su vez en la misma dirección que la corriente exterior. En este caso, se dice que la capa límite está adherida al sólido, puesto que aunque el espesor de esta capa vaya cambiando siempre será tal que δ/c1 (véase la Figura A.2). Sucede sin embargo que, cuando la geometría abre paso al fluido, este se decelera, due/ds<0 , lo que implica la existencia de un gradiente adverso de presión, es decir la presión aumenta aguas abajo conforme el fluido se va decelerando. Aunque en la zona exterior a la capa límite, este efecto se traduce en una disminución en el módulo de la velocidad sin que exista cambio de dirección, en las regiones donde los efectos viscosos son importantes y 50 Capítulo A. Introducción a la Aerodinámica potencial Figura A.2 Distribución de la presión sobre la superficie de un perfil aerodinámico. En aquellas regiones en las que el gradiente de presiones es favorable, la capa límite permance adherida a, siendo su espesor característico δ∼cRe−1/2c . Sin embargo, si el ángulo de ataque es lo suficientemente elevado, el gradiente de presiones por el extradós del perfil puede ser lo suficientemente adverso como para desprender la corriente en esta zona. Esto genera una región de recirculación en la que la presión es prácticamente unifrorme y muy parecida a la del punto de desprendimiento, pd<p∞. donde la velocidad del fluido relativa al cuerpo decrece hasta anularse, este gradiente adverso de presión puede provocar que el fluido en la capa límite pueda cambiar de sentido, comenzando a recircular en contra de la corriente exterior. Si esto ocurre, la capa límite deja de estar adherida a la superficie del cuerpo, y por supuesto deja de tener un espesor muy inferior a las dimensiones características del dominio fluido, formándose una burbuja de recirculación; a este fenómeno se le conoce como desprendimiento de capa límite. Las burbujas de recirculación son características en la zona de sotavento de sólidos romos (cuerpos donde la relación entre la longitud paralela a la corriente incidente y la transversal es de orden unidad) o bien de sólidos fuselados donde el ángulo de ataque α con el que incide la corriente no es pequeño en absoluto. Dado que los valores de la velocidad en esta región son muy inferiores a los existentes fuera de la misma, la presión varía muy poco. Por lo tanto, la presión en el interior de la burbuja es muy parecida a la que tenía el fluido antes de desprenderse, y dado que en el punto de desprendimiento ue>U∞ , la presión en este punto aplicando Euler-Bernoulli será pd−p∞∼−O(ρU2 ∞)<0 . Por otra parte, en la parte de barlovento, donde las partículas se deceleran, se tendrá que la presión es del orden de pb−p∞∼O(ρU2 ∞)>0 , hecho que supone que el coeficiente de resistencia, definido como CD=D/(1/2ρU2 ∞A) , con A el área frontal del objeto, sea de orden unidad. Además, si bien el hecho de que la extensión de las burbujas de recirculación no sea conocido a priori hace que la aerodinámica (hidrodinámica en este caso) de cuerpos romos o fuselados a altos ángulos de ataque no sea abordable mediante técnicas analíticas , cuando el Reynolds es elevado, ni siquiera puede resolverse mediante simulación numérica. Sin embargo, para la mayoría de las aplicaciones se requiere que el coeficiente de resistencia tenga un valor muy inferior a la unidad, lo que nos lleva a eliminar la burbuja de recirculación del problema a través del empleo de cuerpos fuselados a bajos ángulos de ataque. De este modo, la capa límite siempre permanece adherida a la superficie del sólido, por lo que los efectos viscosos pueden ser despreciados en todo el dominio excepto en dicha capa viscosa. Así, las ecuaciones que permiten determinar los campos de presión, velocidad y temperatura pueden escribirse del siguiente modo ∇·v=0,(A.26) ρ∂v ∂t+ρv·∇v=−∇p,(A.27) A.4 Ecuación de la vorticidad y flujos potenciales 51 ρDe Dt=ρcv ∂T ∂t+ρcvv·∇T=−p∇·v⇔ρcv DT Dt=0(A.28) donde en la última ecuación se ha hecho uso de la ecuación de continuidad simplificada para el caso de un líquido incompresible. Nótese que, en este caso, la derivada sustancial de la temperatura es nula, lo que implica que esta se mantiene constante para un observador que se mueve a la velocidad del fluido, v . Este hecho permite desacoplar el problema térmico del mecánico, por lo que nos centraremos en adelante en resolver este último haciendo uso de las Ecuaciones (A.26) - (A.27) . Estas ecuaciones deben ser resueltas junto con las condiciones de contorno a grandes distancias del objeto x→∞,p−p∞→0,v→U∞(A.29) y la única condición de contorno sobre la superficie del objeto, que es la condición de impenetrabilidad: v·ns=0(A.30) Para enunciar esta condición de contorno se considera que la velocidad normal al sólido en las regiones próximas a la capa límite variará muy poco de su valor en la pared, ya que ∆V(yδ)/U∞∼δ/c∼Re−1/2 , por lo que se están cometiendo errores del orden de Re−1/2 al emplear la condición de impenetrabilidad. El mismo orden de error es el que se comete cuando se impone la condición de impenetrabilidad sobre la superficie del cuerpo en lugar de sobre la frontera de la capa límite. A.4 Ecuación de la vorticidad y flujos potenciales Las ecuaciones (A.26) - (A.27) junto con las condiciones de contorno (A.29) - (A.30) permiten determinar el campo de velocidades y de presión de un líquido incompresible que circula alrededor de un cuerpo fuselado a bajos ángulos de ataque. Una simplificación más es posible si se considera que el flujo incidente proviene de una región en la que los valores de la presión y la velocidad son uniformes e iguales a p∞ y U∞ , respectivamente. Al ser el valor de la velocidad del fluido uniforme aguas arriba del objeto, su vorticidad, ω ω ω=∇×v , es idénticamente nula. En efecto, conviene recordar que la vorticidad es una medida de la velocidad angular con la que rotan las partículas adyacentes a un punto material fluido, ω ω ω=2Ω Ω Ω , por lo que en el caso de un flujo uniforme no existe tal rotación de las partículas. Si ahora se hace uso de la identidad vectorial v·∇v=∇v2 2−v×(∇×v) = ∇v2 2−v×ω ω ω(A.31) la ecuación de cantidad de movimiento puede expresarse como ρDv Dt=∂v ∂t+∇v2 2−v×ω ω ω=−∇p(A.32) Si ahora se calcula el rotacional de la ecuación anterior, se tiene en cuenta que ∇·ω ω ω=0 , que ∇×∇ψ=0 siendo ψ cualquier campo escalar, y que los operadores de derivada temporal, ∂/∂t y ∇ pueden conmutar, la ecuación anterior puede simplificarse dando como resultado ∂ω ω ω ∂t+v·∇ω ω ω=Dω ω ω Dt=ω ω ω·∇v−ω ω ω∇·v=ω ω ω·∇v⇐⇒ ⇐⇒ Dω ω ω Dt=ω ω ω·∇v (A.33) donde se ha hecho uso de la ecuación de continuidad, ∇·v=0 . La ecuación anterior expresa que la variación por unidad de tiempo del vector vorticidad siguiendo a una partícula fluida está causado por el término denominado como vortex stretching, ω ω ω·∇v . Este hecho implica que el valor de la vorticidad de una partícula varía a una velocidad que es proporcional al valor instantáneo de la propia vorticidad, por lo que si inicialmente ω ω ω(x→∞) = ∇×U∞ex , se puede concluir que la vorticidad continuará siendo nula en el resto del dominio. Obviamente, el razonamiento seguido parte de que es posible la aplicación de la Ecuación (A.27) , válida a distancias del orden c del objeto donde Rec1 para un líquido incompresible de propiedades constantes. 52 Capítulo A. Introducción a la Aerodinámica potencial En efecto, a distancias del orden del espesor de la capa límite, δ , los efectos de la viscosidad no pueden ser despreciados, y por lo tanto la presencia de vorticidad en estas zonas está asegurada. Teniendo en cuenta, según lo anterior, que el valor de la vorticidad es despreciable en la mayor parte del dominio fluido de nuestro problema, por lo que ∇×v=0 , es posible escribir que v=∇φ , siendo φ una función escalar que en adelante denominaremos potencial de velocidades. Ahora pues, si se introduce este resultado en la Ecuación (A.26), se tiene ∇2φ=0(A.34) por lo la resolución del problema del flujo de un líquido incompresible a altos números de Reynolds alrededor de un objeto fuselado a bajos ángulos de ataque se reduce a hallar la solución de la laplaciana, junto con la ecuación de cantidad de movimiento para el valor de la presión ∇∂φ ∂t+1 2|∇φ|2+p ρ=0=⇒∂φ ∂t+1 2|∇φ|2+p ρ=C(t)(A.35) que no es más que la ecuación de Euler-Bernoulli junto con el término de derivada temporal del potencial (donde se ha hecho uso de la conmutabilidad de los operadores ∇ y ∂/∂t ), y sujeta a las condiciones de contorno x→∞,∇φ→U∞ex,p−p∞→0(A.36) x∈Σs,∇φ·ns=0(A.37) En principio puede parecer que el problema queda perfectamente cerrado de esta forma. No obstante, pese a existir soluciones sencillas de la laplaciana que se expondrán en la próxima sección, se podrá comprobar como estas son infinitas para el problema de flujo potencial alrededor de un cuerpo, quedando indeterminadas hasta que se especifique la condición que permite que el flujo resultante sea físicamente posible. A esta condición se le conoce con el nombre de condición de Kutta-Joukowski y será presentada en la Sección A.7. A.5 Resolución de la ecuación de Laplace mediante superposición de soluciones elementales La linealidad de la ecuación de Laplace en el problema (A.34) - (A.37) permite aplicar el principio de superposición para hallar la solución del problema total como suma de soluciones elementales de la laplaciana, es decir, φ=Σi=N i=1φi(A.38) donde cada φi satisface el problema la Ecuación (A.34) y de forma que el potencial total φ satisface, además, las condiciones de contorno. La mayoría de las superficies sustentadoras, como es el caso de un ala, tienen la peculiaridad de que la longitud transversal principal a la corriente incidente (la envergadura) es mucho mayor que la longitud paralela a la misma (cuerda), por lo que el flujo, en primera aproximación, puede suponerse como bidimensional. Esta simplicación implica que podremos despreciar el término ∂2φ/∂y2 , lo que facilitará la comprensión de los principios de resolución que aquí se presentan. De este modo, a lo largo de esta y las siguientes secciones, se obtendrán soluciones a la laplaciana y se resolverá el problema de un flujo potencial en el caso bidimensional, tomando como sistema de coordenadas el plano x−y con el eje x dirigido según la cuerda del perfil y el eje y perpendicular a este y apuntando hacia el extradós del perfil; nótese que no debe confundirse dicho sistema con los ejes del caso tridimensional, donde y recorre de forma paralela la envergadura del ala. Por otro lado, en el Apéndice B se abordará el problema tridimensional directamente pues así lo requiere la aplicación detallada en el Capítulo 3. A.5.1 Función potencial correspondiente a una corriente uniforme, φ1 Una corriente uniforme que forma un ángulo α con el eje de abcisas x puede expresarse en la forma U∞=U∞cosαex+U∞sinαey , por lo que teniendo en cuenta que ∇xi=∂xi/∂xjej=ei con i=1,2 , se tiene A.5 Resolución de la ecuación de Laplace mediante superposición de soluciones elementales 53 que U∞=U∞cosα∇x+U∞sinα∇y=∇(U∞cosαx+U∞sinαy) = ∇φ1(A.39) por lo que el potencial de velocidades asociado a U∞es φ1=U∞cosαx+U∞sinαy=U∞·x(A.40) A.5.2 Función potencial correspondiente a un manantial o un sumidero, φ2 Considérese la expresión de la ecuación de Laplace en coordenadas polares ∇2φ=1 r ∂ ∂rr∂φ ∂r+1 r2 ∂2φ ∂θ2(A.41) El caso más simple de esta ecuación lo constituye aquel en el que existe independencia de la coordenada azimutal, θ, por lo que φ2=φ2(r)debe satisfacer la ecuación ∇2φ=1 r ∂ ∂rr∂φ ∂r=0→r∂φ2 ∂r=C→φ2=Clnr,(A.42) siendo Cuna constante. Dado que en coordenadas polares ∇φ=∂φ ∂rer+1 r ∂φ ∂θ eθ(A.43) y r=px2+y2,er=xex+yey px2+y2,eθ=−yex+xey px2+y2,(A.44) la expresión del campo de velocidades correspondiente a la solución elemental φ2es v2=∇φ2=er C r=Cxex+yey r2(A.45) Nótese que el valor de la velocidad tiende a infinito cuando r→0 , por lo que existe una singularidad en torno al punto donde se sitúa dicha solución, así como que |v2|→0 cuando r→∞ . A aquellas soluciones donde C>0 se las denomina fuentes mientras que cuando C<0 se denominan sumideros. La razón de esta denominación se halla en que el caudal que circula a través de cualquier superficie Σ que rodee al origen de coordenadas (esto es, a la singularidad) entrará hacia la superficie si C<0 o saldrá de la misma si C>0 . Este hecho permite calcular el valor de la constante como q=ZΣ v·ndσ=Z2π 0 v·nr0dθ=Z2π 0 er·er C r0 r0dθ=CZ2π 0 dθ=2πC=⇒C=q 2π(A.46) Para aquellos casos en los que la fuente o sumidero se halle situada en un punto distinto al origen de coordenadas, por ejemplo x0= (x0,y0) , tan sólo será necesario realizar una traslación en la forma x0=x−x0 , con lo que las soluciones pueden expresarse como φ2=q 2πlnq(x−x0)2+(y−y0)2(A.47) v2=q 2π(x−x0)ex+(y−y0)ey (x−x0)2+(y−y0)2(A.48) Para ilustrar la idea de que la superposición de soluciones elementales también es solución de la ecuación de Laplace, se propone calcular el campo de velocidades generado por una distribución continua de fuentes situadas en el eje horizontal entre x=0 y x=` , siendo el caudal inyectado por la fuente localizada en x0=x0ex igual a dq(x0) = v(x0)dx0 , con v(x0) conocida. El campo de velocidades generado por dicha singularidad será dv(x) = v(x0)dx0 2π(x−x0)ex+yey (x−x0)2+y2(A.49) 54 Capítulo A. Introducción a la Aerodinámica potencial por lo que el campo de velocidades total creado por una superposición contninua de fuentes como la descrita es v(x) = Z` 0 v(x0) 2π(x−x0)ex+yey (x−x0)2+y2dx0(A.50) A.5.3 Función potencial correspondiente a un torbellino, φ3 De forma análoga a la sencillez que spone considerar que una solución de tipo fuente o sumidero es independiente de la coordenada acimutal, es posible considerar el caso contrario de una solución que sólo dependa de dicha coordenada y no de la radial, φ3=φ3(θ). En este caso, la ecuación de Laplace se reduce a ∇2φ3=1 r2 ∂2φ3 ∂θ2=0⇒∂φ3 ∂θ =C⇒φ3=Cθ(A.51) con lo que la velocidad creada por esta solución elemental es v3=∇φ3=∇(Cθ) = er ∂φ3 ∂r+eθ r ∂φ3 ∂θ =eθ C r=C−yex+xey x2+y2(A.52) Esta solución recibe el nombre de torbellino debido a que el campo de velocidades que genera sólo posee componente acimutal. Al igual que en el caso de la fuente, el campo de velocidades creado por un torbellino puntual decrece con el radio, disipándose por completo cuando r→∞ , y existiendo una singularidad en el origen de coordenadas donde se sitúa al torbellino. La constante C puede ser relacionada con el valor de la circulación, Γ , del campo de velocidades a lo largo de un camino C que encierra al eje de coordenadas y que se recorre en sentido antihorario Γ=IC v3·d` ` `(A.53) con lo que tomando un círculo de radio genérico r0 , d` ` `=r0dθeθ , por lo que Γ=R2π 0Cdθ=2πC y de este modo v3=Γ 2πreθ(A.54) φ3=Γ 2πθ(A.55) Finalmente, de forma análoga al caso de las fuentes o sumideros, es posible hallar el campo de velocidades creado por una distribución de torbellinos distribuidos a lo largo del eje x . Llamando dΓ(x0) = v(x0)dx0 a la circulación diferencial creada por un torbellino situado en x=x0 , la expresión de la velocidad creada por esta distribución es análoga a la del caso de la fuente y puede escribirse como v(x) = Z` 0 v(x0) 2π−yex+(x−x0)ey (x−x0)2+y2dx0.(A.56) A.5.4 Función potencial correspondiente a un doblete, φ4 Esta solución elemental surge cuando se pretende calcular el potencial creado por una fuente y un sumidero situados a una distancia 2ε . Ambas singularidades poseen una intensidad idéntica q de tal forma que qε/π→M , en el límite ε→0 , siendo M una constante denominada instensidad del doblete. Para hallar la expresión del potencial correspondiente a un doblete supongamos, sin pérdida de generalidad, que colocamos una fuente de intensidad q en r1=ε,θ=π y un sumidero de intensidad −q en r2=ε,θ=0 . Empleando el princpio de superposición, el potencial de velocidades vendrá dado por π=q 2π[ln||r+εex||−ln||r−εex||](A.57) A.5 Resolución de la ecuación de Laplace mediante superposición de soluciones elementales 55 Puesto que ||r+εex||=p(r+εex)(r+εex) = qr2+2εr·ex+ε2≃ ≃qr2+2εr·ex+Oε2≃r1+2εr·ex r21/2 =r1+2εr·ex r2(A.58) donde se ha tenido en cuenta que, al ser ε1 , √1+ε≃1+1/2ε+O(ε2) y ln(1+ε)≃ε+O(ε2) . Por tanto, ln||r+εex|| =lnr+εr·ex r2(A.59) ln||r−εex|| =lnr−εr·ex r2(A.60) y en el límite ε→0, φ4=Mer r·ex=Mcosθ r(A.61) donde se ha tenido en cuenta que qε/π→M cuando ε→0 . Conocido el potencial es fácil determinar el campo de velocidades creado por el doblete en la forma v4=∇Mer r·ex=−Mer r2er·ex+Meθ r2eθ·ex=M r2(−cosθer−sinθeθ)(A.62) habiéndose tenido en cuenta que ∂er,θ/∂r=0y que ∂er/∂θ =eθ. A.5.5 Potencial creado por la superposición de soluciones elementales Hasta el momento se han hallado 4 soluciones elementales a la Laplaciana en dos dimensiones, lo que nos permite combinarlas linearmente para obtener todo un abanico de soluciones analíticas a diferentes problemas sencillos. Queda fuera del alcance de este apéndice abordar las diferentes soluciones al flujo potencial alrededor de objetos que pueden obtenerse aplicando el principio de superposición con estas soluciones elementales; en [ 19 ] pueden encontrarse multitud de casos resueltos y propuestos que permitan al lector profundizar más en este tema. No obstante, conviene exponer un caso particular que resulta no sólo de gran interés sino que ayudará a comprender el origen de la última condición de contorno al problema (A.34) - (A.37) . Consideremos el flujo resultante de la superposición de una corriente unifrome de módulo U∞ paralela al eje x ( α=0 ) con un doblete y un torbellino, ambos en el origen de coordenadas y con intensidades respectivas U∞a2y−ΓEl potencial resultante de la suma de estas tres soluciones elementales es el siguiente: φ=U∞x+U∞a2cosθ r−Γθ 2π(A.63) cuyo campo de velocidades puede obtenerse a partir de v=∇φ, v=U∞cosθ1−a2 r2er+−U∞sinθ1+a2 r2−Γ 2πreθ(A.64) donde se ha tenido en cuenta que e1=ercosθ−eθsinθ . En la expresión anterior, puede notarse que el valor de la componente radial de la velocidad, vr , es nulo independientemente del valor de θ para la circunferencia de radio r=a centrada en el origen. Por lo tanto, el campo de velociades y el potencial hallados anteriormente satisfancen la ecuación de Laplace sujeta las condiciones de contorno ∇φ→U∞ para r→∞ y n·∇φ=0 sobre el círculo de radio r=a centrado en el origen, con lo que se concluye que dicho campo potencial es solución del problema de flujo potencial de una corriente uniforme que incide sobre un cilindro de radio a. Sin embargo, como puede comprobarse, el campo de la Ecuación (A.63) satisface tanto la ecuación de Laplace como las condiciones de contorno de impenetrabilidad y en el infinito independientemente del valor de la intensidad del torbellino Γ , lo que pone de manifiesto la existencia de infinitas soluciones al problema a 56 Capítulo A. Introducción a la Aerodinámica potencial expensas de proporcionar una condición de contorno que permita reducir el abanico de soluciones a la única físicamente real y posible. A.6 Fuerzas sobre objetos bidimensionales según la teoría potencial Como ya se comentó en la introducción de este apéndice, no es necesario en nuestro caso el cálculo de las fuerzas de sustentación y resistencia sobre el perfil aerodinámico (o el ala, en su caso) ya que para obtener las leyes de escala de los modelos descritos en el Capítulo 2 y Capítulo 3 sólo es necesario conocer el campo de presiones sobre el objeto para poder calcular el gradiente de presión en el punto de interés. Sin embargo, como ha quedado mostrado en el último apartado de la sección anterior, existen infinidad de soluciones al problema (A.34) - (A.37) si no se añade una condición de contorno más. Dado el eminente sentido físico de esta condición, el conocimiento de las fuerzas que actúan sobre el objeto ayudará a obtener lo que se conoce como la condición de Kutta-Joukowski, por ello le dedicaremos un instante a evaluar primero las fuerzas que actúan sobre el cilindro de radio r=a del apartado anterior, extendiendo finalmente dicha expresión para el caso más general de fuerzas sobre cuerpos de geometría arbitraria. A.6.1 Fuerzas sobre un cilindro de radio r=ainmerso en un flujo potencial La resultante de las fuerzas que actúan sobre un cuerpo para el caso del flujo potencial es igual a la resultante de las fuerzas de presión (ya que se despreciaron en su momento los esfuerzos viscosos en la ecuación de cantidad de movimiento). Esta fuerza viene dada por ff s =ZΣs (p−p∞)(−ns)dσ(A.65) siendo ns la normal exterior al sólido. Para hallar la presión en todo punto del dominio fluido, podemos recurrir como ya se sugirió en la sección anterior a la ecuación de cantidad de movmiento (Ecuación (A.35) ), que particularizada para el caso de flujo estacionario resulta en la ecuación de Euler-Bernoulli p+1 2ρv2=p∞+1 2ρU2 ∞(A.66) Si ahora se tiene en cuenta que sobre el cilindro de radio r=ael valor de la velocidad es ∇φ(r=a) = vθeθ=−2U∞sinθ+Γ 4πU∞aeθ(A.67) la presión sobre la superficie del cilindro puede hallarse como p−p∞=1 2ρU2 ∞−v2 θ =1 2ρU2 ∞−2ρU2 ∞sinθ+Γ 4πU∞a2(A.68) la fuerza sobre el cilindro estará dada por la integral ff s =ZΣs (p−p∞)(−ns)dσ=2ρU2 ∞Z2π 0sin2θ+Gsinθcosθex+sinθeyadθ,(A.69) donde G=Γ/(2πU∞a) y donde se ha tenido en cuenta que por el Teorema de Gauss RΣsCnsdσ=0 , con C una constante. Esta ecuación puede reescribirse ahora de la siguiente forma ff s =2aρU2 ∞Z2π 0sin2θcosθ+Gsinθcosθex+sin3θ+Gsin2θeydθ(A.70) A.6 Fuerzas sobre objetos bidimensionales según la teoría potencial 57 y dado que Z2π 0 sin2cosθdθ=0,Z2π 0 sin3θdθ=0,Z2π 0 sinθcosθdθ,y Z2π 0 sin2θdθ=π, (A.71) se tiene ff s =2πaρU2 ∞Gey=ρU∞Γey(A.72) Puesto que el vector resultante de fuerzas que el fluido ejerce sobre el cilindro sólo tiene componenete según la dirección perpendicular a la corriente incidente, se tiene que la teoría potencial predice que el cilindro no ejerce resistencia aerodinámica alguna, lo que se conoce como la paradoja de D 0 Alembert. Evidentemente, este resultado no es físicamente posible si se tienen en cuenta la multitud de experimentos existentes y que pueden consultarse (por ejemplo) en [ 19 ]. La razón de esta aparente contradicción no sólo estriba en la resistencia de fricción causada por la viscosidad que desde un principio se está despreciando, sino además en que en el proceso de resolución del flujo potencial no se contempla en ningún caso la posibilidad de que la capa límite se desprenda, provocando que la presión no se recupere en la zona de sotavento y que por lo tanto el coeficiente de resistencia final sea de orden unidad. A.6.2 Fuerzas sobre objetos bidimensionales de geometría arbitraria inmersos en un flujo potencial El objetivo de esta sección es demostrar que la expresión hallada para la fuerza que el fluido ejerce sobre un cilindro puede extenderse al caso de un objeto de geometría arbitraria, todo ello sin tener que resolver el campo de velocidades potencial. Considerese un volumen de control, Ωc ,como el de la Figura A.3 delimitado por la superficie Σc=Σs∪Σ∞ . La ecuación de cantidad de movimiento aplicada a dicho volumen de control es d dtZΩc ρvdω+ZΣc ρvv ·ndσ=ZΣc (p−p∞)(−nc)dσ= =ZΣf (p−p∞)−nfdσ+ZΣ∞ (p−p∞)(−n∞)dσ= =−ff s +ZΣ∞ (p−p∞)(−n∞)dσ, (A.73) donde se ha tenido en cuenta que ff s =RΣs(p−p∞)(−ns)dσ y que la normal exterior al volumen de control Ωc sobre la superficie Σs es igual y de sentido contrario a la normal exterior sobre el sólido, ns=−nf . En el caso de flujo estacionario y teniendo en cuenta que los flujos de cantidad de movimiento a través de las paredes del sólido son idénticamente nulos, se obtiene que fs=Z∞(p−p∞)(−n∞)dσ−ZΣ∞ ρvv ·n∞dσ(A.74) Antes de proceder al cálculo de las integrales de superficie de la ecuación anterior, se va a calcular el valor de la presión sobre la superficie Σ∞ . En efecto, sobre esta superficie, al estar muy lejos del sólido, la velocidad será prácticamente la velocidad incidente más una perturbación, esto es, x∈Σ∞,v=U∞ex+u0 con u0U∞. La presión en esta región puede determinarse haciendo uso de la ecuación de Euler-Bernoulli, p(x∈Σ∞)−p∞=1 2ρU2 ∞−1 2ρ(U∞ex+u0)·(U∞ex+u0) = =1 2ρU2 ∞−1 2ρU2 ∞+2U∞u0·ex+u0 2= =−ρU∞u0·ex+Ou0/U∞2 (A.75) 64 Capítulo A. Introducción a la Aerodinámica potencial E Et:+ ,, -- ¿\ ,,,-- ,". ,~ \ ll E+ \( ¿ ,_ _,, ' ll E- f Figura A.9 Esquema de las distintas superficies utilizadas para calcular las integrales que involucran a la singularidad. habiéndose tenido en cuenta que l´ ımε→0εlnε→0 . Por otra parte, aplicando el teorema de Gauss y puesto que en el volumen encerrado por la superficie Σs∪Σε−,∇2φ∞=0 y ∇2ψj=0 al no existir singularidades en el interior del mismo, I∞,s=−ZΩs,∞ψj∇2φ∞−φ∞∇2ψjdω=0(A.98) Puede ahora atajarse el cálculo de la integral para el potencial φ en la Ecuación (A.95) sobre la superficie Σe: ZΣeψj∇φ−φ∇ψj·ndσ= l´ ım ε→0Zπ 0 [(lnε)∇φ·nε−(φ/ε)(er)·(−er)]εdθ=πφ xj, (A.99) y, por tanto, la ecuación original puede ser expresada como φx0j,z0j=2cosαx0j+sinαz0j+1 πZΣs−ψj∇φ+φ∇ψ·ndσ+ Γ πZ∞ 0 ∇ψj·(−kcosα+isinα)ds (A.100) donde xj=x0ji+z0jk y se ha tenido en cuenta que la etela sigue la dirección de la corriente incidente. Teniendo en cuenta ahora la condición de impenetrabilidad sobre la superficie del objeto, n·∇φ=0 , y que n·∇ψj=n·x−xj x−xj·x−xj(A.101) la ecuación de Green queda, finalmente, φx0j,z0j=2cosαx0j+sinαz0j+1 πZΣs φ(x,z)x−xj·n x−xj 2dσ+ Γ πZΣex−xj·(−kcosα+isinα) x−xj 2dσ (A.102) donde xj=x0ji+z0jk es el punto sobre la superficie del perfil Σs sobre el que se sitúa el doblete y, por tanto, es un punto fijo a efectos de la integración. Las N ecuaciones de la Ecuación (A.102) (nótese que existe una ecuación para cada posición xj de los N dobletes que se colocarán sobre la superficie del perfil) serán las que permitan resolver numéricamente el flujo sobre un cuerpo de geometría conocida. El primer paso a ejecutar cuando se realiza una implementación numérica consiste en discretizar la geometría. Así, se divide el perfil en N segmentos denominados en adelante paneles. Cada panel k está definido por los nodos k y k+1 y denotaremos a su longitud ∆k . El último nodo A.8 Resolución numérica del flujo alrededor de un perfil aerodinámico: método de los elementos de contorno 65 ( N+1 ) coincide con el nodo 1 , lo que da como resultado los N paneles mencionados. Denotando como rk el vector de posición del nodo k, se puede escribir el punto medio del panel y la longitud del mismo como mk=rk+rk+1 2=xk+xk+1 2i+zk+zk+1 2k=mkxi+mkzk(A.103) ∆k=qrk+1−rk·rk+1−rk=qxk+1−xk2+zk+1−zk2(A.104) Del mismo modo, conocidos los vectores de posición de cada nodo y la longitud de cada panel, pueden encontrarse los vectores tangente y normal al panel kcomo tk=rk+1−rk ∆k (A.105) nk=zk+1−zk ∆k i−xk+1−xk ∆k k(A.106) Conviene definir además la distancia entre dos puntos medios consecutivos, dada por ∆0 k=rm(k+1)x−mkx2+m(k+1)z−mkz2(A.107) Una vez se ha discretizado y definido la geometría, el siguiente paso consiste en colocar el origen singular del doblete j en el punto medio medio del panel j , esto es xj=mj y aproximar el valor del potencial en todo el panel j como una constante igual al valor del potencial en el punto medio de dicho panel. Teniendo en cuenta todo lo anterior, la ecuación integral de contorno puede escribirse como φj=2mjx cosα+mjz sinα+1 π N ∑ k=1 φkZ∆k 0x−mj·nk x−mj 2ds Γ πZ∞ 0x−mj·(−kcosα+isinα) x−mj 2ds (A.108) donde s es la longitud de arco en cada panel, φj=φmj y φk=φ(mk) , lo que proporciona un sistema de N ecuaciones para cada valor de j∈[1,N] . Existen dos integrales que será preciso resolver para formular cada una de las ecuaciones, la primera de ellas, para cada panel k, Z∆k 0x−mj·nk x−mj 2ds(A.109) Si k=j , esto es, elementos diagonales, el valor de la integral es cero, puesto que (x−mk) es perpendicular a nk . En el caso general de k6=j , expresando el vector de posición x sobre el panel k como x=rk+stk , con s∈[0,∆k], la integral queda Z∆k 0xk+stk−mj·nk rk+stk−mj·rk+stk−mjds(A.110) Dado que rk+stk−mj=rk−mj·tktk+rk−mj·nknk+stk(A.111) desarrollando el denominador y realizando la integral se obtiene Z∆k 0xk+stk−mj·nk rk+stk−mj·rk+stk−mjds= =arctan"∆k+rk−mj·tk rk−mj·nk#−arctan"rk−mj·tk rk−mj·nk#(A.112) 66 Capítulo A. Introducción a la Aerodinámica potencial Este resultado puede ser utilizando también para calcular la integral sobre la estela Z∞ 0xbs +s(cosαi+sinαk)−mj·(−kcosα+isinα) xbs +si−mj·xbs +si−mjds= ="π 2sign1−mjxsinα+mjz cosα−arctan 1−mjxcosα−mjz sinα 1−mjxsinα+mjz cosα!# (A.113) donde se ha tenido en cuenta que el vector de posición del borde de salida, que es de donde parte la estela, es xbs =i . Obsérvese que se dispone de N ecuaciones para N+1 incógnitas: los N valores del potencial φk y el valor de la circulación Γ . Se precisa incorporar la condición de Kutta-Joukowski para completar el sistema. Esta condición se implementa simplemente imponiendo que la velocidad en el borde de salida ha de ser igual en módulo y sentido cuando se calcula por el intradós o el extradós, es decir φM−φM−1 ∆0 M−1 =φM+1−φM+2 ∆0 M+1 (A.114) donde se ha tenido en cuenta que el borde de salida está en el nodo M+1 . De este modo, el método de los elementos de contorno o método de Green expuesto hasta este punto, permite hallare el potencial sobbre la superficie del perfil mediante la resolución del siguiente sistema lineal de N+1ecuaciones    ∑N k=1δjk −Rjk−Γfj=2mjx cosα+mz j sinα φM−φM−1 ∆0 M−1−φM+1−φM+2 ∆0 M+1 =0(A.115) para las N+1 incógnitas [φ1,φ2,...,φN,Γ] , que son los valores del potencial en los puntos medios de cada panel y Γ la circulación definida como la integral de línea de la velocidad sobre una curva C que engloba al perfil aerodinámico y está recorrida en sentido horario. Los coeficientes del sistema anteriorse definen como Rjk =     1 π arctan"bs jk +∆k bn jk #−arctan"bs jk bn jk #! k6=j 0k=j      ,(A.116) con bs jk =rk−mj·tkybn jk =rk−mj·nk(A.117) y fj=1 πhπ 2sign 1−mjxsinα+mjz cosα −arctan 1−mjxcosα−mjz sinα 1−mjxsinα+mjz cosα!# (A.118) Así, finalmente, una vez que se dispone del valor del potencial en cada panel, y recordando que la expresión del coeficiente de presión en el caso potencial y utilizando la ecuación estacionaria de Euler-Bernoulli Cp=p−p∞ 1 2ρU2 ∞⇒1−u U∞2 =⇒1−∂φ ∂x2 (A.119) se puede calcular numéricamente el coeficiente de presión como Cp(rk) = 1−φk−φk−1 ∆0 k−12 (A.120) con lo que ya se dispone de todo lo necesario para evaluar el gradiente de presión local en los puntos que se requieran. Apéndice B Aerodinámica de alas de envergadura finita. Vortex-Lattice B.1 Ecuación del potencial y condiciones de contorno En este apéndice se aborda directamente el caso tridimensional de un ala de envergadura, b , finita. Considérese el flujo alrededor de un ala referida al sistema de ejes viento definido en la Figura A.1, con el eje x siguiendo la dirección del flujo U∞, el eje z perpendicular a éste y hacia arriba y el eje y dirigido según el ala derecha del avión. La geometría del ala está dada por las funciones z=z+ w(x,y) y z=z− w(x,y) correspondientes a su extradós e intradós, respectivamente, y quedan determinadas por la forma en planta del ala, la forma de los perfiles de las distintas secciones del ala y el ángulo de ataque que el perfil central forma con la corriente incidente. En el Apéndice A se simplificaron las ecuaciones de Navier-Stokes al realizar varias hipótesis acerca del tipo de movimiento que el fluido describía alrededor del cuerpo. Algunas de las principales hipótesis realizadas son: • Se considera un líquido incompresible de propiedades termodinámicas constantes. Esta hipótesis elimina la densidad como incógnita y la fija como parámetro. Además, simplifica la ecuación de continuidad a ∇·v=0(B.1) • Las fuerzas másicas que actúan son las fuerzas de gravedad, y por lo tanto estas derivan de un potencial, fm=−∇U. • La velocidad característica del líquido, U∞∼O(0.5m/s) , son tales que, junto con el valor característico de la cuerda, c∼O(0.1m) y con los valores de la densidad y viscosidad del agua, ρ≃103kg/m3 y µ=10−3Pa ·s , respectivamente, hacen que el número de Reynolds característico sea tal que Re 1 . Este aspecto implica que se puedan despreciar los efectos viscosos a distancias del orden cdel ala. • A distancias del orden δ∼cRe−1/2 los efectos viscosos son del mismo orden que los efectos de inercia, por lo que en esta región no pueden despreciarse; a esta región se le conoce con el nombre de capa límite viscosa. • Un razonamiento similar puede hacerse en la ecuación de la energía, donde los efectos de conducción son despreciables a distancias del orden c del objeto debido a que el número de Peclet, Pe =RePr con Pr el número de Prandtl definido en el Apéndice A, que surge de comparar el término de calor convectivo frente al conductivo en la ecuación de la energía es tal que Pe 1 . Del mismo modo, a distancias del orden δT∼cPe−1/2 de la superficie del cuerpo, los efectos de conducción han de ser tenidos en cuenta. Además, teniendo en cuenta que se está considerando un líquido incompresible, este resultado permite desacoplar el problema mecánico del térmico. • El flujo incidente es uniforme y por lo tanto carente de vorticidad, o lo que es lo mismo ω ω ω=∇×v=0 , lo que unido a la ecuación de la vorticidad del apéndice anterior da como resultado que el campo de velocidades es irrotacional y deriva de un potencial, v=∇φ. 67 68 Capítulo B. Aerodinámica de alas de envergadura finita. Vortex-Lattice • Este último resultado nos lleva, junto con la ecuación de continuidad para el caso de un líquido incormpresible, a que el problema a resolver es la ecuación de Laplace para el potencial φ ∇2φ=0(B.2) Además de las hipótesis realizadas hasta ahora, y dado que el potencial asociado a la corriente uniforme U∞ es el mismo que en el caso bidimensional, φ∞=U∞·x , con x=xex+yey+zez , se va a resolver la ecuación anterior para el potencial de perturbación, φ0=φ−φ∞ . Este hecho implica que las variables de interés del problema pueden también expresarse como las variables asociadas al potencial infinito más una perturbación, esto es V=U∞ex+uex+vey+wez(B.3) p(x) = p∞+p0(x)(B.4) El hecho de considerar el flujo como la superposición de una solución básica más una perturbación y dado que estas serán pequeñas respecto a las soluciones en cuestión, implica que se pueden linearlizar ciertos términos en la ecuación de Navier-Stokes que permiten escribir de forma más adecuada las condiciones de contorno que se verán a continuación. En efecto, si se considera la ecuación simplificada de cantidad de movimiento para este problema, ρv·∇v=−∇py se linealiza se obtiene ρv·∇v=−∇p→ρU∞ex·(∇v0) = −∇p0→ρU∞ ∂v0 ∂x=−∇p0(B.5) con lo que, proyectando en la dirección exe integrando se obtiene, p0+ρU∞u=C(x,y)(B.6) y dado que en (x,y)→∞ las perturbaciones tienden a cero, se tiene que C(x,y) = 0 , por lo que p0=−ρU∞u . Para este problema, las condiciones de contorno son: •Condición de contorno en el infinito: x→∞,φ→0(B.7) • Sobre la superficie del ala debe cumplirse la condición de impenetrabilidad, es decir V·n± w=0, donde n+ w y n− w denotan las normales unitarias y exteriores al extradós y al intradós, definidos por las ecuaciones z−z+ w(x,y) = 0, y z−z− w(x,y) = 0 respectivamente. Puesto que las normales son proporcionales a ∇(z−z± w) = ez−∂z± w ∂xex−∂z± w ∂yey, la condición de impenetrabilidad puede escribirse como (U∞ex+v0)·∇(z−z± w) = wx,y,z± w−U∞+ux,y,z± w∂z± w ∂x−vx,y,z± w∂z± w ∂y=0(B.8) y teniendo en cuenta que : a) u,vU∞ y b) ∂z± w ∂x∼Oh0 c1 y que ∂z± w ∂y∼Oh0 b1 (ésta última sólo en el caso de que el ála sea plana, puesto que si tuvieramos un ángulo diedro ϒ , entonces ya no se cumpliria exactamente) wx,y,z± w−U∞ ∂z± w ∂x=0 donde, linealizando w teniendo presente que ésta varía en distancias del orden del espesor del ala, tenemos que wx,y,0±=U∞ ∂z± w ∂x(B.9) teniendo que (x,y)εFP,con FP la forma en planta del ala. • Sobre la estela, dado que se trata de una superficie de discontinuidad, debe satisfacerse que la resultante de las fuerzas de presión sobre ella sea nula y la ecuación de continuidad. En efecto, consideremos un cilindro diferencial como volumen de control que engloba a la estela. Aplicando el equilibrio de las B.1 Ecuación del potencial y condiciones de contorno 69 fuerzas de presión sobre todas las superficies y teniendo en cuenta que si el cilindro es diferencial y se hace tender su altura a cero, las superficie lateral no afecta al equilibrio tenemos que p+ E(x,y,z) = p− E(x,y,z)→p+0 E(x,y,z) = p−0 E(x,y,z)(B.10) Recordando ahora que ya habíamos definido para todo el campo fluido una expresión para las perturbaciónes de presión, se tiene: uEx,y,z=0+=u− Ex,y,z=0−∀(x,y)εFP de laestela (B.11) que, a pesar de que es una condición de contorno en velocidad, viene de haber aplicado un equilibrio entre cargas. Esta condición puede recordar a la de Kutta-Yukowsky que se aplica sobre perfiles, es decir, que la velocidad tangencial no varíe en el borde de salida para que no se produzca el rebordeo de la corriente en ese punto). Aplicando ahora la ecuación integral de continuidad al volumen diferencial de control, tenemos que a través de la superficie lateral apenas tendremos flujo por ser un volumen diferencial, mientras que por la superior e inferior ρE+VE+·nE+dΣ+ E=ρE−VE−·nE−dΣ− E(B.12) nE±=±ez(B.13) de donde se obtiene V0 E−·(−ez) = V0 E+·ez(B.14) wx,y,0+=wx,y,0−∀(x,y)εFPde la estela (B.15) El problema por lo tanto puede resumirse, con sus condiciones de contorno                          ∇2·φ=0∀x |∇φ|→0|x|→∞ ∂φ ∂z(x,y,0±) =U∞ ∂z± w ∂xsi(x,y)∈FP ∂φ ∂x(x,y,0+) =∂φ ∂x(x,y,0−) si(x,y)∈E ∂φ ∂z(x,y,0+) =∂φ ∂z(x,y,0−) si(x,y)∈E (B.16) Como ya habíamos adelantado, se trata de un problema de ecuaciones diferenciales en derivadas parciales con las condiciones de contorno apropiadas. Su solución es complicada, pero podemos aprovecharnos (una vez más) de la linealidad del problema y expresar el problema como superposición de otros dos. Recordando que la geometría del ala está dada por las funciones z−z± w(x,y) = 0, para el extradós y el intradós, podemos dividir esta en 2 funciones: una de ellas, zs(x,y), modelará el espesor del ala y la otra, za(x,y) , modelará el ángulo de ataque y la curvatura. De esta forma, en el problema de simétrico o espesor, la resultante de presión será igual en el extradós y en el intradós, y por lo tanto no tendremos sustentación alguna; por contra, en el problema antisimétrico o sustentador, sí tendremos una resultante total del presiones. Veamos ahora como podemos formular ambos problemas. En efecto, el potencial de perturbación φ(x,y,z) admite ser escrito en la forma φ(x,y,z) = φs(x,y,z)+φa(x,y,z)(B.17) siendo φs(x,y,z) = 1 2(φ(x,y,z) +φ(x,y,−z)) (B.18) φa(x,y,z) = 1 2(φ(x,y,z)−φ(x,y,−z)) (B.19) 70 Capítulo B. Aerodinámica de alas de envergadura finita. Vortex-Lattice por lo que puede comprobarse que φs(x,y,z) = φs(x,y,−z)(B.20) φa(x,y,z) = −φa(x,y,−z)(B.21) Si ahora dividimos la geometría de una forma parecida tenemos zs(x,y) = z+ w+z− w 2(B.22) za(x,y) = z+ w−z− w 2(B.23) quedando a su vez que ∂z± w ∂x=∂za ∂x±∂zs ∂x(B.24) De este modo, el problema total se puede dividir en dos subproblemas conocidos como problema espesor y problema sustentador. El primero de ellos es un problema simétrico que puede resolverse mediante una superposición de soluciones elementales de tipo fuente. El segundo es un problema antisimétrico, y es el que considera la línea media cada perfil del ala (tanto si estos tienen curvatura como si no), junto con la corriente U∞ incidiendo con un ángulo de ataque, α . Este es el problema que se requiere resolver para obtener el ángulo de ataque inducido necesario para el cálculo de el gradiente de presión local en el ala rotatoria, por lo que nos centraremos en la resolución de este. Por lo tanto, el problema a resolver es                          ∇2·φa=0∀x |∇φ|→0|x|→∞ ∂φa ∂z(x,y,0±) =U∞ ∂za ∂xsi(x,y)∈FP ∂φa ∂x(x,y,0+) =∂φa ∂x(x,y,0−) =0 si(x,y)∈E ∂φa ∂z(x,y,0+) =∂φa ∂z(x,y,0−) si(x,y)∈E (B.25) donde se ha tenido en cuenta que φax,y,0+=−φx,y,0− en la condición de contorno de velocidad longitudinal en la estela. Una vez identificado pues el problema que se pretende resolver con el fin de hallar el ángulo de ataque efectivo del ala, se va a describir un método numérico conocido como el Método de Vortex-Lattice, ampliamente utilizado hoy día en el cálculo de la sustentación y resistencia inducida de alas durante los procesos de diseño. B.2 Resolución del problema sustentador: método de Vortex Lattice El método de Vortex-Lattice parte de la analogía existente con el problema de flujo potencial y el problema del campo magnético. Es por ello que se describe en primer lugar la analogía del campo de velocidades con la ley de Biot y Savat. Así, se pueden distribuir a lo largo del ala una serie de torbellinos que son soluciones elementales de la Laplaciana como se vio en el Apéndice A y obligarlos a cumplir las condiciones de contorno, obtiendo así un sistema lineal de ecuaciones que se puede resolver numéricamente. B.2.1 Campo de velocidades. Analogía de Biot y Savart La ley de Biot y Savart proporciona la siguiente la expresión del campo magnético alrededor de un hilo infinito por el que circula una intensidad I B=Il µ0I(x0) 4πdl0×x−x0 |x−x0|3(B.26) siendo l , una curva cerrada que rodea al conductor, x e l punto donde queremos medir el campo magnético y x0 la coordenada que reccore el hilo donde situamos la intensidad de corriente. La analogía del campo B.2 Resolución del problema sustentador: método de Vortex Lattice 71 electromagnético con el campo de velocidades aerodinámicas estriba en la definición de la vorticidad como ω=∇×v(B.27) de forma que como el flujo potencial carece de vorticidad, entonces v=∇φ de manera análoga a como sucede con el campo magnético donde el campo eléctrico ( ∇×E=0→E=∇V→∇2V=0 ). Por ello, vamos a proponer como solución al campo de velocidades una expresión del tipo v=∇φ=c(x0)×(x−x0) |x−x0|3=(x−x0)ex+(y−y0)ey+ (z−z0)ez h(x−x0)2+(y−y0)2+ (z−z0)2i3/2(B.28) siendo c(x0)un término constante. Recordando que (x−x0) |x−x0|3=−∇1 |x−x0|(B.29) se tiene v=−c(x0)×∇1 |x−x0|=∇1 |x−x0|×c(x0)(B.30) veamos si el campo propuesto satisface primero nuestra ecuación. En efecto, en el caso incompresible que estamos tratando, la ecuacíón a resolver es la laplaciana, de forma que debe cumplirse que ∇·v=0 . En efecto ∇·v=∇·∇1 |x−x0|×c(x0)=c(x0)·∇×∇1 |x−x0|=0(B.31) donde se ha aplicado la permutabilidad cíclica del producto mixto y que el rotacional de un gradiente es cero por definición. Además sabemos que el campo vdebe satisfacer que ∇×v=0, de modo que ∇×∇1 |x−x0|×c(x0)=∇1 |x−x0|(∇·c(x0))−c(x0)∇21 |x−x0|=0(B.32) donde el primer termino es cero por ser c(x0) un término constante en cada punto x0 , y el segundo también lo es puesto que sabemos de campos electromagnéticos que V=q |x−x0|, y que se cumple que ∇2V=0 . Por lo tanto, al cumplirse las condiciones necesarias para que se satisfaga la ecuación de laplace, ya tenemos un campo vectorial para el campo de velocidades. Para terminar de concretarlo, definiremos el término c(x0) = Γ 4πdl0, por lo que v=Il Γ 4πdl0× x−x0 |x−x0|3!(B.33) Es ahora el momento de ver cómo es el campo de velocidades que crea un hilo de torbellinos, es decir, es el momento de ver cómo resolvemos esta integral. Considerese el segmento definido por PQ, es decir que va desde el punto P a Q. Realizando el producto vectorial del integrando, podemos expresar 4v=Γ 4πdl0|x−x0| |x−x0|3sinβ=Γ 4πdl0 n |x−x0|2sinβ(B.34) siendo β el ángulo que forma (x−x0) con dl0 y n el vector normal a estos dos. Ahora, llamando r al módulo del vector que une los puntos del hilo de torbellinos con el punto donde queremos medir la velocidad, r=|x−x0|,tenemo que la distancia del hilo de torbellinos al punto es d=rsinβ, por lo que tan(π−β) = d l0 =sinβ −cosβ=−tanβ(B.35) por lo que se puede despejar l0como l0=d −tanβ→dl0=d tan2β 1 cos2βdβ=d sin2βdβ(B.36) 72 Capítulo B. Aerodinámica de alas de envergadura finita. Vortex-Lattice pudiendo asi expresar ahora la velocidad como ∆v=Γ 4πnd sin2β sin2β d2sinβdβ=Γ 4πdnsinβdβ(B.37) y, finalmente, integrando entre P y Q, que forman ángulos β1yβ2,tenemos que v=Zβ2 β1 Γ 4πdnsinβdβ=Γ 4πdn(cosβ1−cosβ2)(B.38) si denotamos ahora sP y sQ los vectores que van desde el punto P y Q al punto donde medimos la velocidad x,respectivamente, tenemos que n=sP×sQ sP×sQ xPQ =xP−xQ d=sP×sQ xPQ donde, finalmente aprovechandonos de la propiedad del producto escalar, podemos expresar los cosenos como cosβ1=xPQ ·sP |sP|xPQ cosβ2=xPQ ·sQ sQxPQ de forma que la ecuación quedaría v=Γ 4π sP×sQ sP×sQ 2xPQ· sP |sP|−sQ sQ!(B.39) cuya componente zestá dada por la expresión wPQ (x) = Γ 4πsP×sQz sP×sQ 2xPQ· sP |sP|−sQ sQ!(B.40) que puede escribirse matricialmente como wPQ (x) = ΓgPQ (x)(B.41) siendo de nuevo gPQ (x), un término puramente geométrcio que cuantifica la influencia del anillo de torbellinos PQ sobre el punto x. Una vez que hemos obtenido una expresión para el campo de velocidades, vamos a comentar brevemente aspectos relativos a la implementación numérica. Para ello, dividimos la forma en planta del ala y la estela en los MxN paneles que comentamos al principio de la sección, identificando cada uno por un doble índice (i,j), donde el primer índice indica la posición del panel a lo largo de la cuerda del ala mientras que el segundo lo hace según la envergadura. Cada uno de los paneles los identificaremos con 4 puntos numerados en el sentido de las agujas del reloj. Así, los puntos 2 y 3 de cada panel coinciden con los 1 y 4 del panel que esta justo a la derecha,respectivamente, y los puntos 3 y 4 coinciden con los 2 y 1 del panel posterior (mayor x ). Para los anillos de torbellinos, como recordamos del comienzo de la sección, tendremos torbellinos distribuidos a lo largo del eje x y del eje y , por lo que vamos a situar un anillo cuadrilátero de torbellinos cuyos segmentos laterales conciden con los del panel, el segmento anterior (menor x ) estará situado a un cuarto del panel 4xi j/4 y el segmento posterior coincidirá con el segmento anterior del siguiente panel en la fila. Por otro lado, los puntos de colocación estarán situados a 34xi j/4 , de manera análoga a como hicimos en el caso bidimensional. Los anillos de torbellinos se enumerarán con las letras A,B,C,D, en el mismo sentido que los números del panel. De esta forma las coordenadas de los vértices del cuadrilátero de torbellinos en función de las del panel son: xA=x1+1 4(x4−x1) , xB=x2+1 4(x3−x2), y ahora xC(k,l) = xB(k+1,l) y xD(k,l) = xA(k+1,l) . Dado que la velocidad vertical estaria ahora dada por cuatro hilos de torbellino por cada panel, si queremos calcular la velocidad en un punto de colocación del panel (i,j), sobre él influirán las intensidades de torbellinos de todos los paneles, multiplicadas por su coeficiente de influencia. De modo que B.2 Resolución del problema sustentador: método de Vortex Lattice 73 el sistema de ecuaciones puede escribirse como wi j =U∞ ∂zp(x,y) ∂x(i,j) = M ∑ k=1 N ∑ l=1 Γi jkl gkl (B.42) siendo gkl = (gAB +gBC +gCD +gDA)kl (B.43) los coeficientes de influencia de cada panel (k,l). Para poder implementar el método agruparemos los 4 índices en dos, I y K, siendo I=N∗(i−1)+ j , donde I va desde I=1,...,MxN . y K=N∗(k−1)+l , con K=1,...,MxN. Por otro lado, cabe recordar que vamos a panelar el ala y la estela, pero hasta ahora no hemos hablado nada más que de dividir el ala en MxN paneles. ¿Qué hacemos con la estela? Como ya se ha comentado, la estela es una superficie de discontinuidad, donde, en la zona de extradós, las velocidades tranversales al flujo van desde la punta hasta la raíz y por la zona de intradós lo hacen desde la raíz hasta la punta, representando por tanto una discontinuidad matemática 1 . Si atendemos a la primera de las condiciones de contorno impuestas para la estela en el problema sustentador ux,y,0+=ux,y,−=−ux,y,0+=0∀(x,y)∈E(B.44) ∂φ ∂x(x,y,±)∈E =0→∂ ∂y ∂φ ∂x(x,y,±)∈E!=0→∂ ∂x ∂φ ∂y(x,y,±)∈E!=0(B.45) donde se ha hecho uso del teorema de Schwarz para intercambiar derivadas. La ecuación anterior expresa que las la variación de la velocidad transversal v es nula en la estela, y por lo tanto igual a la discontinuidad que habría en el borde de salida. Por consiguiente, la condición de contorno en la estela vamos a caracterizarla computacionalmente colocando todos los paneles de la fila M de forma que empiecen en la estela y terminen en el infinito aguas abajo. B.2.2 Aspectos relevantes de la implementación numérica Vamos a comentar algunos de los aspectos más de la implementación numérica, como son los parámetros utilizados y la forma de calcular los coeficientes de influencia, entre otros. Se debe comenzar definiendo la geometría, es decir, los puntos xi,i=1 : 4, y los puntos xj,j=A,B,C,D para los hilos de torbellinos. En cuanto a los centros de torbellinos, sabemos que los puntos A y B están a un cuarto de cada panel, pero el centro del torbellino que se situa en el panel ha de estar en el punto medio, es decir, xg=1 2(xA+xB). Algo similar ocurre en los puntos de colocación, xc=1 2x1+3 4(x4−x1)+x2+3 4(x3−x2) . Para la implementación de la ecuación para la velocidad velocidad vertical, debemos definir correctamente los siguientes parámetros: •Vectores del punto de colocación al inicio y final del hilo de torbellino: sPysQy sus módulos sP=xc−xp,yc−yp(B.46) sQ=xc−xQ,yc−yQ(B.47) asQ=qxc−xp2+yc−yp2(B.48) asP=qxc−xQ2+yc−yQ2(B.49) •Dirección de la normal, componente z sP×sQz,y su módulo sP×sQz=xc−xpyc−yQ−xc−xQyc−yp (B.50) 1 Recuérdese que esta discontinuidad es sólo matemática; en la realidad la viscosidad se encarga de hacer que esta velocidad varie de forma continua a través de la estela en un espesor relativo del orden de δE(x0) c0∼qνx0 U∞c2 0 Índice de Tablas 2.1 Dimensiones relevantes del túnel hidrodinámico 16 2.2 Posición de la abcisa, adimensionalizada con la cuerda del perfil, de los orificios de inyeccón de aire. 19 2.3 Rango paramétrico empleado durante la realización de los experimentos. 22 3.1 Medidas del tanque empleado en los experimentos del ala acuática rotatoria 32 3.2 Dimensiones características relevantes del diseño del ala y los orificios de inyección. El valor de g en la última columna representa la distancia de separación entre los diferentes orificios a lo largo de la envergadura 34 3.3 Espacio paramétrico a explorar en la campaña experimental. La última columna representa las posiciones de los orificios que son analizados en cada vídeo, donde rtrepresenta la distancia de cada orificio al eje de rotación, por lo que la velocidad incidente para una velocidad de rotación dada, Ωi, y para un orificio concreto, rtjes Vi j =Ωirtj36 81 Bibliografía [1] J. Rodríguez-Rodríguez, A. Sevilla, C. Martínez-Bazán, and J. M. Gordillo, “Generation of Microbubbles with Applications to Industry and Medicine,” Annual Review of Fluid Mechanics, vol. 47, no. September, pp. 405–429, 2015. [2] K. Ferrara, R. Pollard, and M. Borden, “Ultrasound Microbubble Contrast Agents: Fundamentals and Application to Gene and Drug Delivery,” Annual Review of Biomedical Engineering, vol. 9, pp. 415–447, 2007. [3] A. L. Klibanov, “Microbubble contrast agents: targeted ultrasound imaging and ultrasound-assisted drug-delivery applications.” Investigative radiology, vol. 41, no. 3, pp. 354–362, 2006. [4] M. Postema and A. I. Applications, “Ultrasound-Induced Gas Release from Contrast Agent Microbubbles,” vol. 52, no. 6, pp. 1035–1041, 2005. [5] F. Garcia-Ochoa and E. Gomez, “Bioreactor scale-up and oxygen transfer rate in microbial processes: an overview.” Biotechnology advances, vol. 27, no. 2, pp. 153–76, 2009. [6] D. Rosso, L. E. Larson, and M. K. Stenstrom, “Aeration of large-scale municipal wastewater treatment plants: State of the art,” Water Science and Technology, vol. 57, no. 7, pp. 973–978, 2008. [7] S. Sander, J. Behnisch, and M. Wagner, “Design of fine-bubble aeration systems for municipal WWTPs with high sea salt concentrations,” Water Science and Technology, vol. 75, no. 7, pp. 1555–1563, 2017. [8] K. Terasaka, A. Hirabayashi, T. Nishino, S. Fujioka, and D. Kobayashi, “Development of microbubble aerator for waste water treatment using aerobic activated sludge,” Chemical Engineering Science, vol. 66, no. 14, pp. 3172–3179, 2011. [9] A. Kawahara, M. Sadatomi, F. Matsuyama, and H. Matsuura, “Prediction of micro-bubble dissolution characteristics in water and seawater,” Experimental Thermal and Fluid Science, vol. 33, no. 5, pp. 883–894, 2009. [10] M. Sadatomi, “Performance of a new micro-bubble generator with a spherical body in a flowing water tube,” vol. 29, pp. 615–623, 2005. [11] F. Garcia-Ochoa and E. Gomez, “Prediction of gas-liquid mass transfer coefficient in sparged stirred tank bioreactors,” Biotechnology and Bioengineering, vol. 92, no. 6, pp. 761–772, 2005. [12] J. Eggers and E. Villermaux, “Physics of liquid jets,” Reports on Progress in Physics, vol. 71, no. 3, p. 036601, 2008. [13] A. Sevilla, J. M. Gordillo, and C. Martínez-Bazán, “Bubble formation in a coflowing air–water stream,” Journal of Fluid Mechanics, vol. 530, pp. 181–195, 2005. [14] J. M. Gordillo, a. Sevilla, and C. Martínez-Bazán, “Bubbling in a co-flow at high Reynolds numbers,” Physics of Fluids, vol. 19, no. May, 2007. 83 84 Bibliografía [15] R. Suryo and O. a. Basaran, “Tip streaming from a liquid drop forming from a tube in a co-flowing outer fluid,” Physics of Fluids, vol. 18, no. 2006, 2006. [16] A. M. Gañán-Calvo and J. M. Gordillo, “Perfectly monodisperse microbubbling by capillary flow focusing.” Physical review letters, vol. 87, no. c, p. 274501, 2001. [17] A. Evangelio, F. Campo-Cortés, and J. M. Gordillo, “Pressure gradient induced generation of microbubbles,” Journal of Fluid Mechanics, vol. 778, pp. 653–668, 2015. [18] J. M. Gordillo, A. Sevilla, and F. Campo-Cortés, “Global stability of stretched jets: conditions for the generation of monodisperse micro-emulsions using coflows,” Journal of Fluid Mechanics, vol. 738, pp. 335–357, 2014. [19] J. M. Gordillo and G. Riboux, Introducción a la Aerodinámica Potencial, 2012. Índice alfabético 85