i
Dep. Ingenie ía Ae oespacial y Mecánica de Fluidos
Escuela Técnica Supe io de Ingenie ía
Uni e sidad de Se illa
Se illa, 2019
Desa ollo e implemen ación en Ma lab de
un código numé ico pa a la esolución de
p oblemas de Ae odinámica po encial
linealizada no es aciona ia
Au o : Fe nando Mo eno Pino
Tu o : José Manuel Go dillo A ias de Saa ed a
T abajo Fin de Mas e
Mas e Uni e si a io en Ingenie ía Ae onáu ica
ii
iii
T abajo Fin de Mas e
Mas e Uni e si a io en Ingenie ía Ae onáu ica
Desa ollo e implemen ación en Ma lab
de un código numé ico pa a la esolución
de p oblemas de Ae odinámica po encial
linealizada no es aciona ia
Au o :
Fe nando Mo eno Pino
Tu o :
José Manuel Go dillo A ias de Saa ed a
Dep. de Ingenie ía Ae oespacial y Mecánica de Fluidos
Escuela Técnica Supe io de Ingenie ía
Uni e sidad de Se illa
Se illa, 2019
i
T abajo Fin de Mas e : Desa ollo e implemen ación en Ma lab de un código numé ico pa a la
esolución de p oblemas de Ae odinámica po encial linealizada no es aciona ia
Au o :
Fe nando Mo eno Pino
Tu o :
José Manuel Go dillo A ias de Saa ed a
El ibunal nomb ado pa a juzga el P oyec o a iba indicado, compues o po los siguien es miemb os:
P esiden e:
Vocales:
Sec e a io:
Acue dan o o ga le la cali icación de:
Se illa, 2019
El Sec e a io del T ibunal
i
ii
A mi amilia
A mis maes os
A mis amigos
iii
ix
Ag adecimien os
En p ime luga quie o da las g acias a mi pad e, a mi mad e y a mi he mano, pues han sido la
base en la que apoya me desde mucho an es que es udia á el G ado en Ingenie ía Ae oespacial
y po supues o du an e el Mas e en Ingenie ía Ae onáu ica.
Da g acias ambién a mi compañe a de iaje, Bea, que desde que la conozco ha es ado ahí pa a
apoya me, calma me, y saca me una son isa cada ez que lo necesi aba y sin es ue zo apa en e.
Muchas g acias po odo
Quie o ag adece ambién a odas las pe sonas que me han demos ado que me quie en y sé
que puedo con a con ellas pa a cualquie cosa: Guille, Viole a, Dani, Elena, Ál a o, F an, Alex,
Pablo… la lis a pod ía con a con una g an can idad de nomb es de pe sonas cuya o ma de se
consigue que siemp e piense que odo a a sali bien. G acias po hace me eí y po habe
compa ido conmigo an os momen os. Es un hono o ma pa e del g upo “Ae oma ia” con
odos us edes.
Da las g acias a las pe sonas que me espe aban en el pueblo cada ez que he uel o: mis abuelas
Gonzala y An onia, mi abuelo “Juani”, mi abuelo Paco, que en paz descanse, y al es o de mi
amilia.
A mis amigos, los “emp esa ios” po que aunque es én locos y sean unos payasos, siemp e se án
mis payasos.
Finalmen e, ag adece a mi u o y acompañan e en es e p oyec o, José Manuel Go dillo, pues
con él empezó mi andadu a en el mundo de la ae odinámica y con él cie o mi e apa uni e si a ia,
con un p oyec o que a a dicho campo. G acias po habe es ado ahí odos y cada uno de los
días que he enido dudas y me he ace cado a consul a le y po el a o que ha enido conmigo.
Ha sido un place .
x i
10. Mon aje............................................................................................................................... 110
x ii
ÍNDICE DE FIGURAS
Figu a 1: Geome ía Típica de un Ala en Flecha ..................................................................................... 23
Figu a 2: Supe icies de In eg ación....................................................................................................... 37
Figu a 3: No males en el en o no de la Singula idad ............................................................................. 38
Figu a 4: Disc e ización de la geome ía de Ala y Es ela ......................................................................... 45
Figu a 5: Mallado iangula simé ico en Ala y Es ela ........................................................................... 45
Figu a 6: Cambio de a iable pa a la in eg ación. G á ico explica i o I ................................................... 47
Figu a 7: Ángulos empleados en la In eg ación...................................................................................... 47
Figu a 8: Ángulos empleados en la disc e ización. G á ico explica i o II................................................. 51
Figu a 9: De inición de 𝑟(𝜃) .................................................................................................................. 51
Figu a 10: Mallado en el en o no de la singula idad .............................................................................. 52
Figu a 11: Po encial en el bo de de salida y en la es ela en égimen es aciona io pa a un alo de y dado
............................................................................................................................................................. 56
Figu a 12: Mallado iangula en égimen es aciona io .......................................................................... 56
Figu a 13: Condiciones de con o no pa a el Po encial en égimen es aciona io ..................................... 57
Figu a 14: No ación Ma icial Usada...................................................................................................... 58
Figu a 15: Condición de Ku a pa a el caso no es aciona io. .................................................................. 61
Figu a 16: P opagación de los po enciales en la Es ela no es aciona ia .................................................. 63
Figu a 17: Disc e ización en alas ec angula es ..................................................................................... 65
Figu a 18: Po encial sob e el ex adós de alas ec angula es ................................................................. 66
Figu a 19: Velocidad pe u bada en el eje x sob e el ex adós de alas ec angula es. ............................ 66
Figu a 20: Velocidad pe u bada en el eje x sob e el ex adós de alas ec angula es. Re inamien o ...... 67
Figu a 21: Disc e ización en alas con lecha .......................................................................................... 68
Figu a 22: Po encial sob e el ex adós de alas con lecha ...................................................................... 69
Figu a 23: Velocidad pe u bada en el eje x sob e el ex adós de alas con Flecha. ................................. 69
Figu a 24: E olución de la Pendien e de la Cu a de Sus en ación con el iempo. Ala Rec angula ......... 71
igu a 25: E olución de la pendien e de la cu a de sus en ación con el iempo. ala con lecha ............. 72
Figu a 26: Es ela no Es aciona ia en caso de Ala en Flecha .................................................................... 73
Figu a 27; E olución de la Pendien e de la Cu a de Sus en ación con el Ala gamien o. Ala Rec angula 77
Figu a 28: E olución de la Pendien e de la Cu a de Sus en ación con el Ala gamien o. Ψ=30º ......... 77
Figu a 29: E olución de la Pendien e de la Cu a de Sus en ación con el Ala gamien o.Ψ=45º .......... 78
Figu a 30: Mo imien o del Fluido sob e el ex adós de ae ona es con lecha ........................................ 79
Figu a 31: Velocidad pe u bada en el eje x sob e el ex adós de alas con Flecha Nega i a ................... 79
Figu a 32: Respues a an e un escalón ................................................................................................... 81
Figu a 33: E olución del Coe icien e de Sus en ación en un mo imien o oscila o io .............................. 82
Figu a 34: E olución del coe icien e de sus en ación en un mo imien o oscila o io. (KyP) [2] ............... 82
x iii
xix
ÍNDICE DE TABLAS
Tabla 1; Resul ados pa a Ala Rec angula en Régimen Es aciona io ......................................... 67
Tabla 2: Resul ados pa a Ala Rec angula en Régimen Es aciona io. Mallado Re inado ........... 67
Tabla 3: Resul ados pa a ala con Flecha 30º en égimen es aciona io. .................................... 69
Tabla 4: Resul ados pa a Ala Rec angula en Régimen No Es aciona io ................................... 71
Tabla 5: Resul ados pa a Ala con Flecha 45º en Régimen Es aciona io .................................... 72
Tabla 6: Resul ados ob enidos po dis in os mé odos pa a alas ec angula es con dis in o ala gamien o
............................................................................................................................................... 75
Tabla 7: Resul ados ob enidos po dis in os mé odos pa a alas con Ψ= 30º con dis in o ala gamien o
............................................................................................................................................... 75
Tabla 8: Resul ados ob enidos po dis in os mé odos pa a alas con Ψ= 45º con dis in o ala gamien o
............................................................................................................................................... 75
Tabla 9: E o Come ido en Alas Rec angula es en Régimen Es aciona io ................................ 80
Tabla 10: E o Come ido en Alas con Ψ=30º en Régimen Es aciona io ................................ 80
Tabla 10: E o Come ido en Alas con Ψ=45º en Régimen Es aciona io ................................ 80
xx
21
1. INTRODUCCIÓN
Desde la in ención de la a iación y has a nues os días, la ecnología empleada en es e campo
ha e olucionado conside ablemen e. Si se compa an los a iones de hoy en día con sus ances os
de la segunda gue a mundial o incluso con los an e io es biplanos, se obse a que el a ance es
ab umado . Los campos en los que más se no a es a di e encia a simple is a son la a iónica que
lle an a bo do las ae ona es ac uales, la capacidad que son capaces de anspo a , la
p opulsión, e c…
Sin emba go y aunque a simple is a no se obse e, a odas las ae ona es se les pide 2 cosas
undamen almen e: que uelen y que aguan en las ca gas a las que se en some idas.
La disciplina que se enca ga de la segunda pe ición que se les hace a las ae ona es es la
elas icidad y la esis encia de ma e iales, la cual es udia el compo amien o de la es uc u a de
la ae ona e en e a la acción de las dis in as ue zas a las que se e some ida du an e el uelo
y es ablece las ca gas lími e que es a es capaz de aguan a sin pone en iesgo a la segu idad de
los pasaje os, ipulación y es uc u a.
Po o o lado, la disciplina que se enca ga de que una ae ona e uele es la mecánica de luidos
y más especí icamen e la ae odinámica. Es a disciplina se enca ga de es udia cómo es el
mo imien o del luido al ededo de la ae ona e. Es o pe mi e modela cómo es el campo de
p esiones sob e la ae ona e y po an o conoce las ue zas que iene que sopo a la es uc u a.
Es as dos disciplinas son la base pa a lo que undamen almen e se les pide a las ae ona es. Sin
emba go ac ualmen e se ha conseguido e de o ma an na u al el hecho de que una ae ona e
despega y a e iza, que la p eocupación de las compañías aé eas ecae en inco po a y
desa olla sis emas que consigan una mayo e iciencia del uelo.
En consecuencia, la ae oelas icidad, que es la disciplina que es udia la in e acción an o de la
pa e elás ica como de la ae odinámica, es una disciplina undamen al en el diseño de cualquie
ipo de ae ona es.
Aquí en an en juego dos concep os: e icacia y e iciencia.
La e icacia se de ine cómo la capacidad de log a el e ec o que se desea, mien as que la
e iciencia implica hace lo con el mínimo desap o echamien o de los ecu sos disponibles. Una
ez que la a iación ya se consiguió que uese e icaz, aho a se desa ollan concep os y se in es iga
pa a mejo a su e iciencia. Si se analiza la e iciencia de un uelo de un pun o a o o, se puede
obse a que el uelo más e icien e es el que se desa olla en égimen de c uce o y a se posible
a un alo de coe icien e de sus en ación cons an e (c uise climb).
22
Es e abajo iene como obje i o desa olla un mé odo de cálculo que pe mi a ob ene alo es
eales del coe icien e de sus en ación de alas pa a dis in a geome ía y compa a los con los
esul ados expe imen ales de la li e a u a de o ma que se puedan ex ae conclusiones ace ca
de las en ajas e incon enien es de las di e en es geome ías. Además, en es e abajo no sólo
se a a es udia el mo imien o en égimen es aciona io, sino ambién el mo imien o de la
ae ona e en égimen no es aciona io.
Pa a ello en los dis in os pun os del abajo se an a a a los siguien es emas:
- En p ime luga se a a lle a a cabo la deducción de las ecuaciones del p oblema
ae odinámico y las simpli icaciones e hipó esis necesa ias pa a plan ea la ecuación a
esol e con el mé odo que se a a plan ea .
- A con inuación se p esen a á la disc e ización del p oblema y esolución numé ica. Se
explica án los aspec os más impo an es a ene en cuen a pa a esol e los p oblemas
numé icamen e.
- T as es o, se exponen los esul ados ob enidos pa a alas de dis in a geome ía an o
pa a el caso es aciona io como pa a el caso no es aciona io y se compa a án con algunos
ob enidos de la li e a u a así como o os ob enidos haciendo uso de o o so wa e
come cial.
- Finalmen e se expond án las conclusiones ob enidas del abajo y se p opond án
p oyec os de in es igación u u os que pe mi an mejo a lo y complemen a lo.
23
2. GEOMETRÍA
En es e pun o se an a analiza las ca ac e ís icas de la geome ía sob e la que se a a ealiza
el es udio. Dicha geome ía consis e p incipalmen e en un ala con o ma en plan a ec angula ,
pe o pos e io men e se ealiza la ex ensión a alas con geome ía apezoidal gené ica.
Se conside a á que el ala iene su o ma en plan a con enida en el plano Z = 0.
Las ca ac e ís icas p incipales que de inen a es e ipo de alas son:
- La en e gadu a (b): De inida como la dis ancia en e los bo des ma ginales del ala.
- La cue da en la aíz (c ): De inida como la dis ancia en e el bo de de a aque y el bo de
de salida del ala medida en la aíz de la misma. En nues o caso es la cue da del ala en
el plano Y=0.
- Supe icie ala (S): Supe icie de la o ma en plan a del ala, de inida como la p oyección
del ala sob e el plano Z=0.
- c : es la cue da en las pun as.
- 𝛹: es el ángulo de lecha del ala medido desde el bo de de a aque.
- α: es el ángulo de a aque.
Se ecue da que pa a alas con lecha, la geome ía iene dada po :
𝑏=√AR·𝑆 ; 𝑐𝑟=2 √𝑆 /AR
1+𝐸 ; 𝑐𝑡=𝐸𝑐𝑟
FIGURA 1: GEOMETRÍA TÍPICA DE UN ALA EN FLECHA
24
𝑥𝑎(𝑦)=|𝑦|∙ an(𝛹) 𝑥𝑠(𝑦)=𝑐𝑟+𝑏2 an(Ψ)+𝑐𝑡−𝑐𝑟
𝑏2|𝑦|
Donde se de inen los pa áme os ala gamien o y es echamien o como sigue:
AR=𝑏2
𝑆 𝐸=𝑐𝑟
𝑐𝑡
Y se ha u ilizado la siguien e nomencla u a:
- xa: de ine los pun os x del bo de de a aque del ala.
- xs: de ine los pun os x del bo de de salida del ala.
La Fo ma en plan a del ala se puede desc ibi de o ma gene al con la siguien e unción:
𝐹(𝑥,𝑦,𝑧,𝑡)≡𝑍w(𝑥,𝑦,𝑡)−𝑧=0
El ala es a á compues a po pe iles ae odinámicos esbel os, lo cual se á necesa io pa a la
pos e io simpli icación y linealización de las ecuaciones gene ales de la mecánica de luidos,
con el in de ob ene unas ecuaciones más sencillas y esolubles. La condición de esbel ez de
los pe iles ae odinámicos implica que las a iaciones de la geome ía a lo la go de las
di ecciones longi udinal y ans e sal del ala son pequeñas. Ma emá icamen e es o es:
𝜕𝑍𝑤
𝜕𝑥~∆𝑍𝑐
∆𝑥𝑐~𝑧0
𝑐𝑟≪1
𝜕𝑍𝑤
𝜕𝑦~∆𝑍𝑐
∆𝑦𝑐~𝑧0
𝑏≪1
Siendo z0 un alo de espeso ca ac e ís ico del ala.
De o ma gene al, la ecuación que de ine la geome ía del ala se puede exp esa como la suma
de 2 con ibuciones. 𝑍𝑤(𝑥,𝑦,𝑡)=𝑍𝑐(𝑥,𝑦,𝑡)±𝑍𝑒(𝑥,𝑦)
Donde Zc(x,y, ) ep esen a los pun os de las líneas de cu a u a de cada pe il en cada ins an e
de iempo y Ze(x,y) es la unción que ep esen a el espeso de cada pe il. Que el alo de la
unción de espeso de cada pe il no dependa del iempo, es indica i o de que el pe il no su i á
ninguna de o mación a lo la go del p oblema.
25
3. ECUACIONES Y CONDICIONES DE CONTORNO
3.1. Ecuaciones Gene ales
En es e pun o se plan ea án las ecuaciones que gobie nan el p oblema a pa i de las ecuaciones
gene ales de la Mecánica de Fluidos y se impond án las condiciones de con o no del p oblema.
El mo imien o del ai e al ededo del ala iene egido po las ecuaciones de Na ie -S okes, las
cuales se pueden esc ibi de o ma gene al como se mues a a con inuación.
Ecuación de Con inuidad:
𝜕𝜌
𝜕𝑡+∇·(𝜌𝑣)=0
(3.1.1)
Ecuación de Can idad de Mo imien o:
𝜌𝜕𝑣
𝜕𝑡+𝜌𝑣·∇𝑣=−∇𝑝+∇·𝜏′+𝜌𝑓𝑚
(3.1.2)
Ecuación de la Ene gía:
𝜌𝑐𝑣𝜕𝑇
𝜕𝑡+𝜌𝑐𝑣𝑣·∇𝑇=−𝑝∇·𝑣+𝜏′:∇𝑣+𝑄𝑟+𝑄𝑞+∇·(𝑘∇𝑇)
(3.1.3)
Es e sis ema con o ma un conjun o de 5 ecuaciones con 6 incógni as:
- 𝜌: Densidad del luido
- 𝑣: Velocidad del luido
- P: P esión del luido
- T: Tempe a u a del luido
Pa a ce a el p oblema se debe comple a el sis ema con la ecuación de es ado del luido bajo
es udio, el cual se a a conside a que es un gas pe ec o.
Ecuación de Es ado:
𝑝
𝜌=𝑅𝑔𝑇
(3.1.4)
Donde 𝑅𝑔=𝑅
𝑀𝑚=𝑐𝑝−𝑐𝑣=𝑐𝑣(𝛾−1). Siendo 𝑅=8,314 𝐽
𝑚𝑜𝑙·𝐾 la cons an e uni e sal de los
gases y 𝑀𝑚 la masa mola del gas.
Con es o, se iene un sis ema no lineal de ecuaciones en de i adas pa ciales con 6 ecuaciones y
6 incógni as al que hab á que p opo ciona le las condiciones de con o no e iniciales adecuadas
pa a cada p oblema.
32
ℎ=𝐶𝑝𝑇=𝐶𝑝𝑝
𝑅𝑔𝜌=(𝐶𝑝
𝐶𝑣)𝑝
(𝑅𝑔
𝐶𝑣)𝜌=𝛾
𝛾−1𝑝
𝜌=𝑎2
𝛾−1
Se iene que:
𝜕𝜙
𝜕𝑡+|∇𝜙|2
2+𝑎2
𝛾−1=𝐶(𝑡)
En es e p oyec o se conside a á únicamen e el lími e 𝑀∞≪1 lo que pe mi e simpli ica las
ecuaciones aún más. En es e lími e 1
𝜌Δ𝜌<< 1. Lo que implica que el luido puede se conside ado
como incomp esible. En es e caso la ecuación de la con inuidad se simpli ica a:
∇·𝑣′
=0
Los é minos siguien e pueden se desp eciados. Véase Go dillo y Riboux, 2012 [1]
∂ρ
∂ ; 𝑣′
·∇𝜌≪𝜌∇·𝑣′
Pa a el caso 𝜌≈𝑐𝑡𝑒, la ecuación escala pa a el po encial queda:
𝜕𝜙
𝜕𝑡+|∇𝜙|2
2+𝑝
𝜌=𝐶(𝑡)
(3.4.2)
Se obse a que la ecuación de can idad de mo imien o ha pasado a se una ecuación escala
pa a el po encial de elocidades. 𝐶(𝑡) es una cons an e que es únicamen e una unción del
iempo. Conocido el po encial de elocidades, se puede u iliza la ecuación an e io pa a
de e mina la dis ibución de p esiones.
En el caso de 𝜌≈𝑐𝑡𝑒, la sus i ución de la elocidad en unción del po encial en la ecuación de
con inuidad p opo ciona que el po encial de elocidades sa is ace la ecuación de Laplace.
∇2𝜙=0
Con la hipó esis de lujo incomp esible se es a ían come iendo e o es de 𝑂(𝑀2), donde
𝑀=𝑣/𝑎 deno a el núme o de Mach. Dicha demos ación se puede encon a en [1]. Es a
hipó esis es álida pa a alas olando a un Mach<0.3.
Debido a que la ecuación de Laplace es lineal, se puede aplica el mé odo de supe posición de
soluciones elemen ales pa a halla la solución del p oblema gene al.
El p oblema gene al queda de inido po :
∇2𝜙=0
|𝑥|→∞, 𝜙→𝑈∞𝑥
(3.4.3)
33
|𝑥|∈Σ𝑠, 𝑛
·(∇𝜙−𝑣𝑝
)=0
𝐶𝑜𝑛𝑑𝑖𝑐𝑖𝑜𝑛 𝑑𝑒 𝐾𝑢𝑡𝑡𝑎−𝐽𝑜𝑢𝑘𝑜𝑤𝑠𝑘𝑖
Donde 𝑣𝑝
es la elocidad e ical del pe il en caso de encon a nos en égimen no es aciona io
y la condición de Ku a-Joukowski es la que pe mi e de e mina la solución eal a nues o
p oblema de en e las in ini as soluciones exis en es que iene el p oblema del lujo po encial
al ededo de un obje o.
3.5. Ecuaciones Linealizadas
A con inuación se a a p ocede a linealiza las ecuaciones a esol e . Se a a conside a que el
ala es á o mada po pe iles muy esbel os y se mue e a ángulos de a aque pequeños. Es o
pe mi e supone que la capa lími e no se a a desp ende en ningún pun o del ala y en
consecuencia, dicho ala pe u ba poco la co ien e inciden e. Es o quie e deci , que las
pe u baciones en las a iables del luido, o iginadas po la p esencia del ala, son de un o den
de magni ud in e io a su alo aguas a iba.
Las hipó esis a conside a son:
𝛼≪1
ℎ0≪𝑐
𝜕𝑧𝑝
𝜕𝑥≪1, 𝜕𝑧𝑝
𝜕𝑦≪1
(3.5.1)
Donde ℎ0 es el espeso ca ac e ís ico del pe il. Es as hipó esis jun o con la de núme o de
Reynolds muy ele ado, nos pe mi en asegu a que la capa lími e no se desp ende y que la
pe u bación en las a iables luidas es pequeña.
Con es as hipó esis los campos de elocidades, p esión, densidad y empe a u a de las a iables
luidas se puede eesc ibi como la suma de su alo aguas a iba más una pe u bación
o iginada po la p esencia del ala (3.5.2):
𝑣(𝑥)=𝑈∞
+𝑣′(𝑥), |𝑣′(𝑥)|≪𝑈∞
𝑝(𝑥)=𝑝∞+𝑝′(𝑥), 𝑝′(𝑥)≪𝑝∞
(3.5.2)
Donde 𝑣′(𝑥), 𝑝′(𝑥), 𝜌′(𝑥), 𝑇′(𝑥) son los campos pe u bados de elocidad, p esión, densidad y
empe a u a.
El po encial de elocidades se puede descompone igualmen e como la suma del po encial de
elocidades en el in ini o más el po encial de elocidades pe u badas (3.5.3), que se á la
incógni a de nues o p oblema.
𝜙(𝑥)=𝜙∞+𝜙′(𝑥)
(3.5.3)
34
Debido a la linealidad de la ecuación de Laplace, el po encial de elocidades pe u badas debe
de cumpli (3.5.4):
∇2𝜙′=0
(3.5.4)
La ecuación de can idad de mo imien o en é minos del po encial de elocidades pe u badas
debe de cumpli (3.5.5)
𝜌𝜕𝜙′
𝜕𝑡+𝜌𝑈∞𝜕𝜙′
𝜕𝑥+𝑝′=0
(3.5.5)
La condición de impene abilidad puede simpli ica se g acias a la linealización de las ecuaciones:
𝑛
·(∇𝜙−𝑣𝑝
)=∇𝐹𝑒,𝑖
|∇𝐹𝑒,𝑖|·(∇𝜙−𝑣𝑝
)=0→∇𝐹𝑒,𝑖·(∇𝜙−𝑣𝑝
)
Donde 𝐹𝑒,𝑖 es la unción de los pun os del ex adós y del in adós del ala. Se puede eesc ibi
pues la ecuación de la siguien e o ma:
(𝑒3
−𝜕𝑧𝑝
𝜕𝑥𝑒1
−𝜕𝑧𝑝
𝜕𝑦𝑒2
)·((𝜕𝜙′
𝜕𝑧−𝜕𝑧𝑝
𝜕𝑡)𝑒3
+(𝑈∞+𝜕𝜙′
𝜕𝑥)𝑒1
+𝜕𝜙′
𝜕𝑦𝑒2
)=0
Desp eciando los é minos de segundo o den, se llega a (3.5.6):
𝜕𝜙′
𝜕𝑧(𝑥,𝑦,𝑧=0)=𝑤′(𝑥,𝑦,𝑧=0)=𝑈∞𝜕𝑧𝑝
𝜕𝑥(𝑥,𝑦,𝑡)+𝜕𝑧𝑝
𝜕𝑡(𝑥,𝑦,𝑡)
(3.5.6)
Es a condición se debe ía de impone es ic amen e sob e el ala, pe o debido a que los pe iles
son esbel os y los ángulos que amos a maneja son pequeños, dicha condición se impone sob e
la o ma en plan a del ala (FP). G acias a es as hipó esis, es a ap oximación no supone un e o
impo an e. Véase Go dillo y Riboux, 2012 [1].
El p oblema linealizado iene pues desc i o po las siguien es ecuaciones:
∇2𝜙′=0
|𝑥|→∞, 𝜙′→0
|𝑥|∈𝐹𝑃,𝜕𝜙′
𝜕𝑧(𝑥,𝑦,𝑡)=𝑈∞𝜕𝑧𝑝
𝜕𝑥(𝑥,𝑦,𝑡)+𝜕𝑧𝑝
𝜕𝑡(𝑥,𝑦,𝑡)
𝐶𝑜𝑛𝑑𝑖𝑐𝑖𝑜𝑛 𝑑𝑒 𝐾𝑢𝑡𝑡𝑎−𝐽𝑜𝑢𝑘𝑜𝑤𝑠𝑘𝑖
(3.5.7)
3.6. Condición de Ku a-Joukowky
En es a apa ado se explica la condición de Ku a-Joukowsky y sus implicaciones en el p oblema
ae odinámico.
35
Dicha condición es ablece que de las in ini as soluciones que iene el mo imien o de un luido
al ededo del ala, es aquella pa a la cual el lujo no ebo dea el bo de de salida. Es a condición
ija el alo de la ci culación (Γ), pe mi iendo elegi de en e odas las soluciones aquella que se
da en la ealidad.
La condición de con o no que se debe de impone en el bo de de salida se ob iene de es a
condición, jun o con la condición de que la di e encia de p esiones en e es ados e in adós
debe de se nula a lo la go del bo de de salida po que no exis e sólido que sopo ase dicha
di e encia en el caso de que la hubiese.
En el ex adós y en el in adós del bo de de salida se debe de cumpli que:
𝜌𝜕𝜙+
𝜕𝑡 +𝜌𝑈∞𝜕𝜙+
𝜕𝑥 +𝑝+=0
𝜌𝜕𝜙−
𝜕𝑡 +𝜌𝑈∞𝜕𝜙−
𝜕𝑥 +𝑝−=0
Si se es an ambas ecuaciones y se ienen cuen a que la p esión en el ex adós y en el in adós
son iguales: 𝜕(𝜙+−𝜙−)
𝜕𝑡 +𝑈∞𝜕(𝜙+−𝜙−)
𝜕𝑥 =0
Teniendo en cuen a la an isime ía del campo de elocidades en el p oblema sus en ado ,
Go dilllo y Riboux, 2012, [1], 𝜙+=−𝜙−, se puede eesc ibi la condición an e io pa a
cualquie a de las supe icies ex adós o in adós. Sup imiendo el supe índice + en el po encial,
se iene inalmen e.
𝜕𝜙
𝜕𝑡+𝑈∞𝜕𝜙
𝜕𝑥=𝐷𝜙
𝐷𝑡=0
(3.6.1)
Es a ecuación es la que se a a impone en el bo de de salida. Es o quie e deci que la a iación
po unidad de iempo del po encial de elocidades siguiendo a la pa ícula luida es nula.
En el caso es aciona io, es a condición se aduce en que en el bo de de salida hay un pun o de
emanso (en el caso de bo de de salida anguloso), o bien que la elocidad en el bo de de salida
iene la misma di ección y sen ido an o en el ex adós como en el in adós (caso de bo de de
salida en e oceso).
En el caso es aciona io la condición an e io queda:
𝜕𝜙
𝜕𝑥=0
(3.6.2)
36
37
4. SOLUCIÓN GENERAL DEL PROBLEMA
Se de ine el po encial 𝜓0 de una solución elemen al de la ecuación de Laplace pa a el caso
idimensional co espondien e a una uen e localizada en el pun o (𝑥0,𝑦0,𝑧0) como:
𝜓0=1
√(𝑥−𝑥0)2+(𝑦−𝑦0)2+(𝑧−𝑧0)2
(4.1)
Se puede comp oba que es e ipo de solución cumple la ecuación de Laplace ∇2𝜓0=0.
El po encial de elocidades pe u badas, debe de cumpli ambién la ecuación de Laplace
∇2𝜙′=0.
Reco demos que se iene que cumpli que:
∇2𝜙′=0, 𝑣′=∇𝜙′, → ∇2𝑣′=0
Donde se ha hecho uso de la in e cambiabilidad de las de i adas.
Debido a es o se debe de cumpli la siguien e ecuación:
𝜓0∇2𝑣′−𝑣′ ∇2𝜓0=0
(4.2)
A con inuación se p ocede a in eg al es a ecuación en un olumen de con ol Ω𝑐 limi ado po la
supe icie Σ𝑐=Σ∞∪Σ𝜖∪Σ𝑎𝑙𝑎∪Σ𝑒𝑠𝑡𝑒𝑙𝑎+∪Σ𝑒𝑠𝑡𝑒𝑙𝑎−
. Dichas supe icies se obse an en 2D en
la siguien e igu a:
∫(𝜓0∇2𝑣−𝑣 ∇2𝜓0)𝑑𝜔=
Ω𝑐∫∇·(𝜓0∇2𝑣−𝑣 ∇2𝜓0)𝑑𝜔=
Ω𝑐
FIGURA 2: SUPERFICIES DE INTEGRACIÓN
38
=∫𝜓0𝑛
·∇𝑣′−𝑣′𝑛
·∇𝜓0 𝑑𝜎=0
Σ𝑐
La in eg al a lo la go de Σ∞ es nula, y las in eg ales a lo la go de Σ𝜖 y Σ𝑎𝑙𝑎 pueden esc ibi se de
la siguien e o ma:
∫𝜓0𝑛
·∇𝑣′−𝑣′𝑛
·∇𝜓0𝑑𝜎=
Σ𝜖∪Σ𝑎𝑙𝑎
=∫𝜓0𝑛
·∇𝑣′−𝑣′𝑛
·∇𝜓0 𝑑𝜎+
Σ𝑎𝑙𝑎 ∫𝜓0𝑛
𝑒−·∇𝑣′−𝑣′𝑛
𝑒−·∇𝜓0 𝑑𝜎+
Σ𝜖−
+∫𝜓0𝑛
·∇𝑣′−𝑣′𝑛
·∇𝜓0 𝑑𝜎+∫𝜓0(−𝑛
𝑒−)·∇𝑣′−𝑣′(−𝑛
𝑒−)·∇𝜓0𝑑𝜎=
Σ𝜖−
Σ𝜖
=𝐼𝑎𝑙𝑎+𝐼𝜖
Donde Σ𝜖 es la semies e a supe io que odea la singula idad y cuya no mal posi i a se de ine
hacia el in e io de la es e a, y Σ𝜖− es la semies e a in e io que odea la singula idad y cuya
no mal posi i a se de ine hacia el ex e io de dicha es e a. La in eg al 𝐼𝜖 se puede calcula de la
siguien e mane a.
𝐼𝜖=lim
ϵ→0∫(1ϵ(−𝑒𝑟)·∇𝑣′(𝑥)−𝑣′(𝑥)(𝑒𝑟)·(−𝑒𝑟)
𝜖2)𝜖2𝑑Ω=−4𝜋𝑣′(𝑥0
)
4𝜋
0
(4.3)
Donde en es e caso 𝑑Ω=𝑑S/ 2 es el di e encial de ángulo sólido.
T as es e cálculo la ecuación de G een puede esc ibi se de la siguien e o ma:
4𝜋𝑣′(𝑥0)=∫𝜓0𝑛
·∇𝑣′−𝑣′𝑛
·∇𝜓0𝑑𝜎
Σ𝑎𝑙𝑎∪Σ𝑒𝑠𝑡𝑒𝑙𝑎+∪Σ𝑒𝑠𝑡𝑒𝑙𝑎−
(4.4)
P oyec ando sob e el ec o uni a io 𝑘
.
4𝜋𝑤
′(𝑥0
)=∫𝜓0𝑛
·∇𝑤
′−𝑤
′𝑛
·∇𝜓0𝑑𝜎
Σ𝑎𝑙𝑎∪Σ𝑒𝑠𝑡𝑒𝑙𝑎+∪Σ𝑒𝑠𝑡𝑒𝑙𝑎−
(4.5)
Debido a que la an isime ía 𝑤′(𝑥,𝑦,𝑧=0+)=𝑤′(𝑥,𝑦,𝑧=0−) y a que las no males al ala y a
la es ela son ap oximadamen e 𝑘
y -𝑘
pa a el ex adós y pa a el in adós, se iene pues que:
∫𝑤′(𝑥,𝑦,𝑧=0 )𝑛
·∇𝜓0𝑑𝜎=0
Σ𝑒𝑠𝑡𝑒𝑙𝑎+∪Σ𝑒𝑠𝑡𝑒𝑙𝑎−
FIGURA 3: NORMALES EN EL ENTORNO DE LA SINGULARIDAD
39
∫𝑤′(𝑥,𝑦,𝑧=0 )𝑛
·∇𝜓0𝑑𝜎=0
Σ𝑎𝑙𝑎+∪Σ𝑎𝑙𝑎−
Luego se llega a la siguien e ecuación
4𝜋𝑤
′(𝑥0)=∫𝜓0𝑛
·∇𝑤
′𝑑𝜎=
Σ𝑎𝑙𝑎∪Σ𝑒𝑠𝑡𝑒𝑙𝑎+∪Σ𝑒𝑠𝑡𝑒𝑙𝑎−
=−∫𝜓0𝜕𝑤′+
𝜕𝑧 𝑑𝜎++∫𝜓0𝜕𝑤′−
𝜕𝑧 𝑑𝜎−
Σ𝑎𝑙𝑎−∪Σ𝑒𝑠𝑡𝑒𝑙𝑎−
Σ𝑎𝑙𝑎+∪Σ𝑒𝑠𝑡𝑒𝑙𝑎+
Aplicando la ecuación de con inuidad an o en el ex adós como en el in adós se iene:
∇·𝑣′+
=𝜕𝑢′+
𝜕𝑥 +𝜕𝑣′+
𝜕𝑦 +𝜕𝑤′+
𝜕𝑧 =0 → 𝜕𝑤′+
𝜕𝑧 =−𝜕𝑢′+
𝜕𝑥 −𝜕𝑣′+
𝜕𝑦
∇·𝑣′−
=𝜕𝑢′−
𝜕𝑥 +𝜕𝑣′−
𝜕𝑦 +𝜕𝑤′−
𝜕𝑧 =0 → 𝜕𝑤′−
𝜕𝑧 =−𝜕𝑢′−
𝜕𝑥 −𝜕𝑣′−
𝜕𝑦
Dado que la placa pa e del eposo, la ci culación Γ es idén icamen e nula en cualquie supe icie
ce ada que englobe al ala como a la es ela según el Teo ema de Bje knes-Kel in, el cual dice:
𝐷Γ
𝐷𝑡=0
(4.6)
Tomando como supe icie, una compues a po el ala y la es ela, se iene:
Γ=∫ ∫ 𝛾(𝑥,𝑦,𝑡) 𝑑𝑥 𝑑𝑦=∫ ∫ (𝑢′+−𝑢′−) 𝑑𝑥 𝑑𝑦=∫ ∫ 𝜕(𝜙+−𝜙−)
𝜕𝑥 𝑑𝑥 𝑑𝑦
𝑥𝑒𝑠𝑡
0
𝑏
2
−𝑏
2
𝑥𝑒𝑠𝑡
0
𝑏
2
−𝑏
2
𝑥𝑒𝑠𝑡
0
𝑏
2
−𝑏
2
Γ=∫(𝜙+(𝑥𝑒𝑠𝑡,𝑦,𝑡)−𝜙−
𝑏
2
−𝑏
2(𝑥𝑒𝑠𝑡,𝑦,𝑡))𝑑𝑦=0
Donde el subíndice “es ” hace e e encia a la es ela y 𝑥𝑒𝑠𝑡 es el inal de la es ela y la a iable
𝛾(𝑥,𝑦,𝑡), es la densidad de ci culación. Se ha u ilizado que el po encial es 0 en el bo de de
a aque, condición que se explica á pos e io men e, cuando se analicen las condiciones de
con o no del p oblema.
Como es a in eg al es nula pa a odo ins an e de iempo, de ello se deduce que el in eg ando
no puede depende del iempo. Es o implica que la a iación con el iempo del po encial en la
coo denada donde acaba la es ela es nula, y po an o las elocidades ho izon ales son iguales.
U ilizando las siguien es igualdades ob enidas de la de i ación de un p oduc o:
𝜓0𝜕𝑢′
𝜕𝑥=𝜕(𝜓0𝑢′)
𝜕𝑥 −𝑢′𝜕𝜓0
𝜕𝑥
40
𝜓0𝜕𝑣′
𝜕𝑦=𝜕(𝜓0𝑣′)
𝜕𝑦 −𝑣′𝜕𝜓0
𝜕𝑦
Y sus i uyendo en la ecuación de G een, se iene que:
4𝜋𝑤
′(𝑥0)=∫(𝜕(𝜓0𝑢′+)
𝜕𝑥 −𝑢′+𝜕𝜓0
𝜕𝑥 +𝜕(𝜓0𝑣′+)
𝜕𝑦 −𝑣′+𝜕𝜓0
𝜕𝑦)𝑑𝜎+−
Σ𝑎𝑙𝑎+∪Σ𝑒𝑠𝑡𝑒𝑙𝑎+
−∫(𝜕(𝜓0𝑢′−)
𝜕𝑥 −𝑢′−𝜕𝜓0
𝜕𝑥+𝜕(𝜓0𝑣′−)
𝜕𝑦 −𝑣′−𝜕𝜓0
𝜕𝑦)𝑑𝜎−
Σ𝑎𝑙𝑎−∪Σ𝑒𝑠𝑡𝑒𝑙𝑎−
Pa a es e azonamien o se a a ene en cuen a que las in eg ales de supe icie se pueden hace
indis in amen e p ime o en la di ección x y luego en la di ección y o ice e sa. Tal como se ha
esc i o la exp esión an e io se obse a que los é minos son iguales en ambas in eg ales,
simplemen e cambiando la no ación de ex adós a in adós. En onces:
- La suma de las in eg ales del p ime é mino es 0, debido a que al hace la in eg al en
di ección x se iene que en el bo de de a aque dichas elocidades son 0, y a que al inal
de la es ela se cumple 𝑢′+=𝑢′− .
- La suma de las in eg ales del e ce é mino es 0 debido a que al hace la in eg al en
di ección y se debe de ene en cuen a que las elocidades 𝑣′ son an isimé icas
espec o al plano z=0, siendo en onces dichas elocidades idén icas pe o de signo
con a io en los bo des ma ginales.
Eliminando es as in eg ales de la exp esión an e io y haciendo uso de nue o de la an isime ía
del campo de elocidades se llega a:
2𝜋𝑤
′(𝑥0)=−∫(𝑢′+𝜕𝜓0
𝜕𝑥+𝑣′+𝜕𝜓0
𝜕𝑦)𝑑𝜎+
Σ𝑎𝑙𝑎+∪Σ𝑒𝑠𝑡𝑒𝑙𝑎+
Ecuación que inalmen e queda:
𝑤
′(𝑥0)=𝑈∞𝜕𝑧𝑝
𝜕𝑥(𝑥0,𝑦0,𝑡)+𝜕𝑧𝑝
𝜕𝑡(𝑥0,𝑦0,𝑡)=−1
2𝜋∫𝑣′(𝑥)·∇𝜓0 𝑑𝜎+
Σ𝑎𝑙𝑎+∪Σ𝑒𝑠𝑡𝑒𝑙𝑎+
(4.7)
Es a es la solución gene al del p oblema ae odinámico idimensional. Se obse a que dicha
exp esión se co esponde con la ecuación de Bio y Sa a . La a iable 𝑣′(𝑥) indica las
componen es de la elocidad en el plano Z=0. El lec o in e esado puede encon a más
in o mación en Go dillo y Riboux, 2012, [1].
De es a o ma se obse a que la componen e e ical del campo de elocidades pe u bado en
cualquie pun o del dominio luido puede se exp esada como la supe posición del campo de
elocidades c eado po una dis ibución con inua de o bellinos si uados sob e la supe icie del
ala y de la es ela, y o ien ados según la en e gadu a, con in ensidad po unidad de longi ud
41
2𝑢′+(𝑥,𝑦), y una dis ibución con inua de o bellinos si uados sob e ala y es ela, de in ensidad
po unidad de longi ud 2𝑣′+(𝑥,𝑦) o ien ados en la di ección del eje x.
T abajos an e io es a an sob e la esolución de es a ecuación in eg al. Sin emba go, dichos
es udios demues an que el sis ema que se ob iene es ines able desde el pun o de is a
numé ico. Véase Flo es Caballe o, 2016 [4]. El mé odo que se a a implemen a co ige es a
ines abilidad median e la ans o mación de dicha ecuación in eg al, ya que en el in eg ando
apa ece á como incógni a el alo del po encial pe u bado en luga de su g adien e.
48
De es a o ma el denominado de la ecuación se puede desa olla como sigue:
|𝑅
−𝑅0
| =|(𝑅
−𝑅1
)+(𝑅1
−𝑅0
)| =|𝑟𝑒𝑟
+(𝑅1
−𝑅0
)| =
=|𝑟𝑒𝑟
+(𝑅1
−𝑅0
)| =(𝑟2+|𝑅1
−𝑅0
|2+2𝑟𝑒𝑟
·(𝑅1
−𝑅0
))1/2=
=((𝑟+𝑒𝑟
·(𝑅1
−𝑅0
))2+|𝑅1
−𝑅0
|2−(𝑒𝑟
·(𝑅1
−𝑅0
))2)1/2=
=(𝑈2+𝑎2)1/2
Donde x y a ienen dadas po las siguien es exp esiones (5.8) y (5.9):
𝑈=𝑟+𝑒𝑟
·(𝑅1
−𝑅0
)
(5.8)
𝑎2=|𝑅1
−𝑅0
|2−(𝑒𝑟
·(𝑅1
−𝑅0
))2
(5.9)
Es a nomencla u a se á u ilizada más adelan e.
U ilizando las nue as a iables 𝜃 y , se llega a la siguien e conclusión.
𝑥−𝑥1=𝑟𝑐𝑜𝑠(𝜃)
𝑦−𝑦1=𝑟𝑠𝑖𝑛(𝜃)
Sus i uyendo en la exp esión an e io , desa ollando el denominado de la in eg al y u ilizando
la de inición de las a iables U y 𝑎2 se llega a lo siguien e.
𝑤𝑖=∫𝜙(𝑥𝑖)
|𝑅𝑖
−𝑅0
|3𝑑𝜎𝑖
Σi= 𝜙1∫[ 1
(𝑈2+𝑎2)3/2+𝐶1𝑟𝑐𝑜𝑠(𝜃)
(𝑈2+𝑎2)3/2+𝐷1𝑟𝑠𝑖𝑛(𝜃)
(𝑈2+𝑎2)3/2]𝑑𝜎𝑖
Σi+
+ 𝜙2∫[ 𝐶2𝑟𝑐𝑜𝑠(𝜃)
(𝑈2+𝑎2)3/2+𝐷2𝑟𝑠𝑖𝑛(𝜃)
(𝑈2+𝑎2)3/2]𝑑𝜎𝑖
Σi+
+ 𝜙3∫[ 𝐶3𝑟𝑐𝑜𝑠(𝜃)
(𝑈2+𝑎2)3/2+𝐷3𝑟𝑠𝑖𝑛(𝜃)
(𝑈2+𝑎2)3/2]𝑑𝜎𝑖
Σi
Donde U=U( ), y los lími es de la in eg al an desde =0, has a = (𝜃). Desa ollando la in eg al
doble y sacando los é minos que no dependen de ue a de la in eg al se llega a lo siguien e.
𝜙1∫ [∫ 𝑟𝑑𝑟
(𝑈2+𝑎2)3/2
𝑟(𝜃)
0+(𝐶1𝑐𝑜𝑠(𝜃)+𝐷1𝑠𝑖𝑛(𝜃))∫𝑟2𝑑𝑟
(𝑈2+𝑎2)3/2
𝑟(𝜃)
0]𝑑𝜃+
𝜃3
θ2
49
+ 𝜙2∫ [(𝐶2𝑐𝑜𝑠(𝜃)+𝐷2sin (𝜃))∫𝑟2𝑑𝑟
(𝑈2+𝑎2)3
2
𝑟(𝜃)
0]𝑑𝜃+
𝜃3
θ2
+ 𝜙3∫ [(𝐶3𝑐𝑜𝑠(𝜃)+𝐷3𝑠𝑖𝑛(𝜃))∫𝑟2𝑑𝑟
(𝑈2+𝑎2)3/2
𝑟(𝜃)
0]𝑑𝜃
𝜃3
θ2
Como se puede obse a , las in eg ales a ealiza son simila es en e ellas. Se a a pasa a
analiza las una a una indi idualmen e.
La p ime a se mues a a con inuación:
𝑆1=∫𝑟𝑑𝑟
(𝑈2+𝑎2)3/2
𝑟(𝜃)
0=∫𝑈−𝑒𝑟
·(𝑅0
−𝑅1
)
(𝑈2+𝑎2)3/2 𝑑𝑈
𝑈(𝜃)
𝑈0
Donde se ha hecho el cambio de a iable mos ado an e io men e. Ope ando con es a exp esión
se ob iene lo siguien e:
∫𝑈−𝑒𝑟
·(𝑅0
−𝑅1
)
(𝑈2+𝑎2)3/2 𝑑𝑈
𝑈(𝜃)
𝑈0=∫𝑈
(𝑈2+𝑎2)3/2𝑑𝑈
𝑈(𝜃)
𝑈0−𝑒𝑟
·(𝑅0
−𝑅1
)∫ 𝑑𝑈
(𝑈2+𝑎2)3/2
𝑈(𝜃)
𝑈0
𝑆1=−𝑆0−𝑒𝑟
·(𝑅0
−𝑅1
)𝑆3
Donde:
𝑆0=−∫𝑈
(𝑈2+𝑎2)3
2𝑑𝑈
𝑈(𝜃)
𝑈0=1
((𝑈(𝜃))2+𝑎2)1
2−1
(𝑈02+𝑎2)1
2
𝑆3=∫𝑑𝑈
(𝑈2+𝑎2)3/2
𝑈(𝜃)
𝑈0=1
𝑎2[𝑈(𝜃)
(𝑈(𝜃)2+𝑎2)1
2−𝑈0
(𝑈02+𝑎2)1
2]
A con inuación se a a analiza el o o ipo de exp esión que apa ece en la in eg al.
𝑆2=∫𝑟2𝑑𝑟
(𝑈2+𝑎2)3/2
𝑟(𝜃)
0=∫(𝑈−𝑒𝑟
·(𝑅0
−𝑅1
))2
(𝑈2+𝑎2)3/2 𝑑𝑈=
𝑈(𝜃)
𝑈0
=∫𝑈2−2𝑈𝑒𝑟
·(𝑅0
−𝑅1
)+(𝑒𝑟
·(𝑅0
−𝑅1
))2
(𝑈2+𝑎2)3/2 𝑑𝑈=
𝑈(𝜃)
𝑈0
=∫𝑈2−2𝑈𝑒𝑟
·(𝑅0
−𝑅1
)+(𝑒𝑟
·(𝑅0
−𝑅1
))2±𝑎2
(𝑈2+𝑎2)3/2 𝑑𝑈=
𝑈(𝜃)
𝑈0
Donde en el úl imo paso se ha sumado y es ado el e mino de 𝑎2. A con inuación, desa ollando
la in eg al y u ilizando las in eg ales que se han calculado an e io men e, se ob iene.
50
𝑆2=∫𝑑𝑈
(𝑈2+𝑎2)1/2+2𝑈𝑒𝑟
·(𝑅0
−𝑅1
)𝑆0+
𝑈(𝜃)
𝑈0[2(𝑒𝑟
·(𝑅0
−𝑅1
))2−|𝑅0
−𝑅1
|2]𝑆3
Po comodidad y al igual que se ha hecho an e io men e, se a a u iliza la siguien e
nomencla u a.
𝑆4=∫𝑑𝑈
(𝑈2+𝑎2)1/2=𝐿𝑛(√𝑎2+𝑈(𝜃)2+𝑈(𝜃)
√𝑎2+𝑈02+𝑈0)
𝑈(𝜃)
𝑈0
𝑆2=𝑆4+2𝑈𝑒𝑟
·(𝑅0
−𝑅1
)𝑆0+[2(𝑒𝑟
·(𝑅0
−𝑅1
))2−|𝑅0
−𝑅1
|2]𝑆3
Sus i uyendo las de iniciones de 𝑆1 y 𝑆2 en la ecuación in eg al pa a el iángulo bajo es udio se
end ía que:
𝑤𝑖= 𝜙1∫[𝑆1+(𝐶1𝑐𝑜𝑠(𝜃)+𝐷1sin (𝜃))𝑆2]𝑑𝜃+
𝜃3
θ2
+ 𝜙2∫[(𝐶2𝑐𝑜𝑠(𝜃)+𝐷2sin (𝜃))𝑆2]𝑑𝜃+
𝜃3
θ2
+ 𝜙3∫[(𝐶3𝑐𝑜𝑠(𝜃)+𝐷3sin (𝜃))𝑆2]𝑑𝜃
𝜃3
θ2
A con inuación es impo an e eco da que 𝑈=𝑈(𝑟(𝜃)), donde 𝑈(𝑟) es conocida, pe o aún
hay que de e mina 𝑟(𝜃). Pa a ello hay que u iliza el siguien e iángulo.
51
Donde se puede e que 𝑟(𝜃) es la dis ancia desde un pun o del lado 3-2 al é ice 1 del
iángulo. Aplicando elaciones en e los dis in os ángulos, se pueden ob ene los ángulos
in e io es de la egión indicada.
Aplicando el Teo ema del seno y despejando se ob iene la exp esión 𝑟(𝜃) buscada.
𝑟(𝜃)
sin (𝜋+𝜃2−𝜃4)=|𝑅2
−𝑅1
|
sin(𝜃4−𝜃)
𝑟(𝜃)=sin (𝜃4−𝜃2)|𝑅2
−𝑅1
|
sin(𝜃4−𝜃)
FIGURA 8: ÁNGULOS EMPLEADOS EN LA DISCRETIZACIÓN. GRÁFICO
EXPLICATIVO II
FIGURA 9: DEFINICIÓN DE 𝒓(𝜽)
52
Una ez se haya hecho la in eg al en la a iable 𝑟, se p ocede a la in eg ación en la a iable 𝜃, la
cual se hace de o ma numé ica median e un mé odo de colocación.
Una ez hecha la in eg ación pa a la p ime a ecuación, ya se ienen los alo es de los
coe icien es de in luencia que acompañan a los po enciales de elocidades. Es os po enciales
an colocados en los é ices de los iángulos en los que se ha di idido el ala y la es ela.
En el caso en que se enga que in eg a en el en o no de la singula idad, se ealiza á el siguien e
mallado:
Y se supond á que el po encial es de la siguien e o ma:
𝜙′(𝑥)=𝜙1′=𝜙0′
Es deci , en el en o no de la singula idad (8 iángulos) el po encial es cons an e e igual al
po encial de la singula idad. Es o e i a la eno me ines abilidad numé ica que end ían las
ecuaciones si se dejase e oluciona el po encial como en el es o de los casos. La in eg al no es
di e gen e en el en o no de la singula idad sólo si el po encial es cons an e.
La condición de con o no de impene abilidad disc e izada queda de la siguien e o ma:
𝑤(𝑥𝑖
)=[𝐴1,1,𝐴1,2,…𝐴1,𝑖,…𝐴1,(𝑁+2)(𝑀+2)]·
[
𝜙1,1
𝜙1,2
…
𝜙1,𝑖
…
𝜙1,(𝑁+2)(𝑀+2)
]
+
+[𝐸1,1,𝐸1,2,…𝐸1,𝑖,…𝐸1,(𝑁+2)(𝑀 𝑒𝑠𝑡𝑒𝑙𝑎)]·
[
𝜙𝐸1,1
𝜙𝐸1,2
…
𝜙𝐸1,𝑖
…
𝜙𝐸1,(𝑁+2)(𝑀 𝑒𝑠𝑡𝑒𝑙𝑎)
]
FIGURA 10: MALLADO EN EL ENTORNO DE LA
SINGULARIDAD
53
Donde 𝜙𝑖,𝑗 indica el po encial en el nodo i,j y 𝐴𝑖,𝑗 y 𝐸𝑖,𝑗 son los co espondien es coe icien es de
in lencia, calculados empleando las exp esiones analí icas de las in eg ales halladas
an e io men e.
54
55
6. PROBLEMA ESTACIONARIO. RESOLUCIÓN Y
CONDICIONES DE CONTORNO
Una ez ealizada la o mulación de la condición de impene abilidad en su o ma gene al, en
es e apa ado se a a pasa a expone las condiciones de con o no pa a el p oblema es aciona io
y su esolución.
Pa a ce a el p oblema, se deben de expone un o al de ecuaciones igual al núme o de
incógni as, con el obje i o de que la ma iz de coe icien es sea cuad ada e in e ible.
Cómo se ha is o an e io men e que en cualquie pun o del dominio luido donde no haya
sólido, se debe cumpli que el po encial de elocidades pe u badas debe de se 0 al se lo en el
in ini o aguas a iba, en pa icula , sob e el bo de de a aque y sob e los bo des ma ginales se
debe de cumpli que, en i ud de la ecuación de Eule -Be nouilli linealizada, que exp esa que la
de i ada sus ancial del po encial pe u bado debe de se 0, se concluye que dicho po encial es
0 en los bo des de a aque y ma ginales.
Bo de de a aque y bo des ma ginales:
𝜙′(𝑥𝑏𝑎,𝑦)=𝜙′(𝑥,𝑏2)=𝜙′(𝑥,−𝑏2)=0
Po o a pa e, en el bo de de salida, la condición de Ku a-Joukowski exp esa en el caso
es aciona io que:
Bo de de salida:
𝜕𝜙′
𝜕𝑥=0→𝜙(𝑥𝑏𝑠−1,𝑦)=𝜙(𝑥𝑏𝑠,𝑦)=𝜙(𝑥𝑏𝑠+1=𝑥𝐸1,𝑦)
Es a condición implica que pa a una línea de coo denada y, el alo del po encial de elocidades
en el bo de de salida en esa línea y es el mismo en el bo de de salida, en la disc e ización en x
an e io al bo de de salida y a lo la go de oda la es ela. G á icamen e se puede e en la igu a
siguien e:
56
Se obse a que en como el po encial es el mismo an o en el bo de de salida, como en el pun o
an e io y como en oda la es ela a lo la go de cada línea y, no iene sen ido ealiza un mallado
a lo la go de la es ela pa a el caso es aciona io. Luego en es e caso, el mallado es mucho más
simple:
FIGURA 11: POTENCIAL EN EL BORDE DE SALIDA Y EN LA ESTELA EN RÉGIMEN ESTACIONARIO PARA UN VALOR
DE Y DADO
FIGURA 12: MALLADO TRIANGULAR EN RÉGIMEN ESTACIONARIO
57
A la ho a de mon a las ecuaciones es e ecin o se puede subdi idi de la siguien e o ma
Únicamen e se debe de impone la condición de con o no de impene abilidad en la zona
in e io al ecuad o ama illo, es deci , un o al de N*M ecuaciones. El po encial en el bo de de
salida y en la es ela ( ecuad o azul) debe de se el mismo que en los pun os an e io es al bo de
de salida, y el po encial en los bo des de a aque y ma ginales ( ecuad os e des) es 0.
Así pues se iene un o al de N*(M+1) incógni as sob e el ala, y N*(M+1) ecuaciones di ididas
en N*M condiciones de impene abilidad, más N condiciones de con o no de Ku a sob e el
bo de de salida. Es o pe mi e o ma una ma iz cuad ada e in e ible con la que esol e el
p oblema.
La ecuación queda ía de la siguien e o ma gene al pa a el p oblema es aciona io (se omi en los
é minos que acompañan a po enciales nulos po la condición de con o no).
𝑤′(𝑥𝑖
)=𝑤𝑖′=[𝐴𝑖,2,𝐴𝑖,3,…𝐴𝑖,𝑀+2,…𝐴𝑖,(𝑁+1)(𝑀+2)]·
[
𝜙2,2
𝜙2,3
…
𝜙2,𝑀+2
…
𝜙2,(𝑁+1)(𝑀+2)
]
+
+[𝐸𝑖,2,𝐸𝑖,3,…𝐸𝑖,𝑖,…𝐸𝑖,(𝑁+1)·]·
[
𝜙𝐸2
𝜙𝐸2
…
𝜙𝐸𝑖
…
𝜙𝐸(𝑁+1)
]
FIGURA 13: CONDICIONES DE CONTORNO PARA EL POTENCIAL EN RÉGIMEN ESTACIONARIO
64
ins an e de iempo simplemen e median e la “con ección” aguas debajo de los po enciales
ob enidos en el bo de de salida.
En el anexo se p opo ciona el p og ama de Ma lab que esuel e el p oblema es aciona io y no
es aciona io pa a geome ías gené icas.
65
8. SOLUCIÓN EN RÉGIMEN ESTACIONARIO.
8.1. Alas Rec angula es
A con inuación se a a expone la solución pa a alas ec angula es en égimen es aciona io, y
pa a dis in os ala gamien os. Pa a ob ene una solución asequible es necesa io alcanza un
equilib io en e el núme o de di isiones del mallado y el iempo compu acional, de mane a que
se in en e alcanza la mejo solución (o una solución con p ecisión su icien e) en un iempo de
cálculo azonable.
Pa a alcanza es e obje i o es p eciso eco da que la dis ibución de p esiones en el ex adós
del ala su e unos g adien es muy acen uados en las p oximidades del bo de de a aque y de los
bo des ma ginales. Es e enómeno lle a a pensa que es necesa io e ina el mallado en es as
zonas, con obje i o de ob ene unos alo es azonables del coe icien e de sus en ación y de la
pendien e de la cu a de sus en ación (Se supond á en odo momen o a iación lineal del
coe icien e de sus en ación con el ángulo de a aque). El mallado en el eje x e y se ob iene con
las siguien es exp esiones ( eniendo en cuen a el núme o de di isiones):
𝑚𝑦=0.5𝑏·cos((𝑗−1)𝜋
𝑁+1) 𝑐𝑜𝑛 𝑗∈[1,𝑁+2]
𝑚𝑥=(1−cos((𝑖−1)𝜋
2(𝑀−1))) 𝑐𝑜𝑛 𝑖∈[1,𝑀]
El mallado pa a el caso de alas en lecha es simila .
Se mues a un ejemplo de cómo se ían las di isiones del mallado. Se ha p esen ado el ejemplo
pa a que se obse e el ala gamien o del ala y se pueda in ui la longi ud de la es ela (la es ela
FIGURA 17: DISCRETIZACIÓN EN ALAS RECTANGULARES
66
llega ía más lejos, pe o la idea con es a imagen es mos a al lec o las dimensiones del
p oblema).
En la igu a an e io se obse a cómo el mallado es más ino en los pun os del bo de de a aque
y ma ginales pa a cap a bien los g adien es.
Resol iendo el p oblema pa a un ángulo de a aque de 5º, un ala gamien o de 5 (b=5, c=1) y
es echamien o nulo, se ob ienen las siguien es dis ibuciones del po encial sob e la supe icie
del ala.
Si se ealiza la de i ación del po encial espec o a la coo denada x, se puede ob ene el campo
de elocidades sob e el ex adós.
Es as dos igu as mues an el eno me g adien e de elocidades que exis e an o en los bo des
ma ginales, como en el bo de de a aque y po an o jus i ican el mallado ealizado.
La pendien e de la cu a de sus en ación ob enida pa a es a con igu ación es de:
FIGURA 18: POTENCIAL SOBRE EL EXTRADÓS DE ALAS RECTANGULARES
FIGURA 19: VELOCIDAD PERTURBADA EN EL EJE X SOBRE EL EXTRADÓS DE ALAS RECTANGULARES.
67
𝐶𝐿𝛼
Nx
Ny
AR
3,8741
15
21
5
TABLA 1: RESULTADOS PARA ALA RECTANGULAR EN RÉGIMEN ESTACIONARIO
A con inuación se ha ealizado el mismo cálculo pe o aumen ando el mallado has a un alo de
Nx=25 y Ny=41 pa a obse a las di e encias.
Se mues a a con inuación una igu a del campo de elocidades pa a es e caso, de o ma que se
ap ecie mejo el g adien e que se p oduce en el bo de de a aque y en los bo des ma ginales.
Los esul ados pa a es e caso son:
𝐶𝐿𝛼
Nx
Ny
AR
3,9096
25
41
5
TABLA 2: RESULTADOS PARA ALA RECTANGULAR EN RÉGIMEN ESTACIONARIO. MALLADO REFINADO
El iempo compu acional ha sido mayo pa a es e caso, y la mejo a se no a en la e ce a ci a
signi ica i a.
Se puede obse a que la pendien e de la cu a de sus en ación iende a aumen a con o me el
mallado se hace más ino. Es o se debe a que se cap a mejo las e oluciones de las a iables
luidas en el en o no de los bo des ma ginales y del bo de de a aque, pe mi iendo que la
in eg ación de luga a un esul ado más p eciso.
De aquí en adelan e se ealiza á el mallado con los siguien es alo es de Nx=15 y Ny=8*b+1. La
azón de es os alo es es po que cuando se esuel a el p oblema no es aciona io, el iempo
FIGURA 20: VELOCIDAD PERTURBADA EN EL EJE X SOBRE EL EXTRADÓS DE ALAS RECTANGULARES.
REFINAMIENTO
68
compu acional aumen a eno memen e pa a mallados muy g andes si se quie en ob ene
alo es que mejo an poco en p ecisión con espec o a es os. De o ma que con es e mé odo
aco amos el e o a la e ce a ci a signi ica i a mien as man enemos un iempo azonable.
La azón de la ecuación Ny=8*b+1, es iba en que en el análisis de alas con g an ala gamien o,
se debe án hace más di isiones en la coo denada y debido a que al aumen a la en e gadu a,
ambién lo hacen los paneles en los que se di ide el ala en esa di ección y po an o se pe de ía
p ecisión si no se aumen a el mallado.
8.2. Alas con Flecha
En es e apa ado se a a la solución del p oblema es aciona io pa a el caso de alas con lecha.
De es a o ma, dando alo es de S, AR, 𝛹 y E, queda de inida el ala.
Pa a el caso de alas en lecha, se p obó a u iliza un mallado simila al caso de alas ec angula es,
sin emba go la p ecisión adop ada en compa ación con los alo es ípicos de la li e a u a no e a
su icien e.
La azón es iba en que las alas en lecha poseen un g adien e de elocidades en las
p oximidades de la cue da aíz del ala, que esul aba se la zona donde el mallado an e io es
menos p eciso (elemen os más g andes, cap aban mal la e olución).
Pa a la comp obación de que es e código e a e ec i o, puede comp oba el lec o in e esado
que si se u iliza el código pa a el ala en lecha con las pa icula idades de E=1 y ángulo de lecha
nulo, se uel e a ob ene los alo es que se ob u ie on pa a el ala con o ma en plan a
ec angula .
Pa a sol en a es e p oblema se ecu ió a un e inado ambién po dicha zona simila a lo que
se había hecho pa a los bo des de a aque y ma ginales. Se mues a el mallado pa icula
ealizado pa a las alas en lecha:
FIGURA 21: DISCRETIZACIÓN EN ALAS CON FLECHA
69
El po encial de elocidades sob e es e ipo de alas iene la siguien e o ma:
De i ando el po encial, se puede ob ene el campo de elocidades pe u badas sob e la
supe icie de ex adós, el cual se mues a en la siguien e igu a:
En es a igu a se puede obse a el g adien e de elocidades que se p oduce en las p oximidades
de la cue da aíz del ala.
𝐶𝐿𝛼
Nx
Ny
AR
3,5368
15
41
5
TABLA 3: RESULTADOS PARA ALA CON FLECHA 30º EN RÉGIMEN ESTACIONARIO.
Se obse a que el alo de la pendien e de la cu a de sus en ación pa a un ala en lecha es
meno que pa a el ala ec angula , como e a de espe a .
FIGURA 22: POTENCIAL SOBRE EL EXTRADÓS DE ALAS CON FLECHA
FIGURA 23: VELOCIDAD PERTURBADA EN EL EJE X SOBRE EL EXTRADÓS DE ALAS CON FLECHA.
70
71
9. SOLUCIÓN EN RÉGIMEN NO ESTACIONARIO.
9.1. Alas Rec angula es
A con inuación se a a mos a los esul ados ob enidos pa a alas ec angula es en égimen no
es aciona io. Pa a ello, se a a ealiza el mismo análisis, con el mismo mallado en el ala, pe o
pa a égimen no es aciona io.
El inc emen o de iempo Δ𝑡 debe de se in e io a 𝑐/𝑈∞, pa a que el e ec o no es aciona io se
obse e co ec amen e.
El iempo de cálculo se ha es ingido a 200 𝑐/𝑈∞, con un Δ𝑡 de 0.25𝑐/𝑈∞, en e cada
i e ación. Las ecuaciones han sido in eg adas en el iempo u ilizando un mé odo de Eule de
p ime o den.
El alo ob enido se mues a a con inuación jun o con la cu a que indica la e olución que ha
seguido el alo de la pendien e de la cu a de sus en ación desde el ins an e inicial has a el
inal pa a el caso de un ala ec angula con elación de aspec o AR=5 y suponiendo que la
elocidad adimensional en el in ini o es una unción escalón en el iempo.
𝐶𝐿𝛼
Nx
Ny
AR
3,8788
15
41
5
TABLA 4: RESULTADOS PARA ALA RECTANGULAR EN RÉGIMEN NO ESTACIONARIO
Cabe esal a el hecho de que la pendien e de la cu a de sus en ación alcanza el alo ob enido
en égimen es aciona io an e io men e.
FIGURA 24: EVOLUCIÓN DE LA PENDIENTE DE LA CURVA DE SUSTENTACIÓN CON EL TIEMPO. ALA
RECTANGULAR
72
9.2. Alas con Flecha
Pa a las alas en lecha se ha seguido un p ocedimien o simila al de las alas ec angula es, con
la pa icula idad de hace un mallado más ino en la zona ce cana a la cue da aíz del ala.
Una ez se iene implemen ado el código pa a la esolución no es aciona ia de alas
ec angula es, simplemen e se modi ica la geome ía pa a ealiza el es udio del
compo amien o no es aciona io pa a alas en lecha.
Se ha ealizado de nue o una es icción en el iempo de cálculo de 200 𝑐/𝑈∞, con un Δ𝑡 de
0.25𝑐/𝑈∞. Es o da paneles de longi ud 0.25 en la es ela y da alo es lo su icien emen e p ecisos
aunque con un iempo compu acional ele ado.
𝐶𝐿𝛼
Nx
Ny
AR
3,5303
15
41
5
TABLA 5: RESULTADOS PARA ALA CON FLECHA 45º EN RÉGIMEN ESTACIONARIO
Pa a que el lec o enga una idea de las dimensiones que iene la es ela no es aciona ia analizada
pa a el case de un ala en lecha se adjun a la igu a siguien e:
FIGURA 25: EVOLUCIÓN DE LA PENDIENTE DE LA CURVA DE SUSTENTACIÓN CON EL TIEMPO. ALA CON FLECHA
73
Cada di isión en el eje x co esponde a un ins an e de iempo. En o al se ienen 200 ins an es
de iempo en di isiones de 0.25 s.
FIGURA 26: ESTELA NO ESTACIONARIA EN CASO DE ALA EN FLECHA
80
El e o ob enido en e el cálculo es aciona io y el ob enido del g á ico se mues an a
con inuación:
ALA RECTANGULAR
AR
𝐶𝐿𝛼𝐸𝑠𝑡𝑎𝑐𝑖𝑜𝑛𝑎𝑟𝑖𝑜
KyP
E o (%)
4
3,5656
3,6
0,95
5
3,8741
3,9
0,66
6
4,1321
4,2
1,61
7
4,3236
4,3
0,54
TABLA 9: ERROR COMETIDO EN ALAS RECTANGULARES EN RÉGIMEN ESTACIONARIO
ALA FLECHA 30º
AR
𝐶𝐿𝛼𝐸𝑠𝑡𝑎𝑐𝑖𝑜𝑛𝑎𝑟𝑖𝑜
KyP
E o (%)
4
3,2911
3,35
1,75
5
3,5368
3,7
4,41
6
3,7222
3,85
3,31
7
3,8712
4
3,22
TABLA 10: ERROR COMETIDO EN ALAS CON 𝚿=𝟑𝟎º EN RÉGIMEN ESTACIONARIO
ALA FLECHA 45º
AR
𝐶𝐿𝛼𝐸𝑠𝑡𝑎𝑐𝑖𝑜𝑛𝑎𝑟𝑖𝑜
KyP
E o (%)
4
2,922
3
2,6
5
3,0942
3,2
3,30625
6
3,2197
3,35
3,8895
7
3,3231
3,5
5,054
TABLA 11: ERROR COMETIDO EN ALAS CON 𝚿=𝟒𝟓º EN RÉGIMEN ESTACIONARIO
Po lo que se puede obse a que el mé odo desa ollado iene un e o máximo del 5% pa a el
caso de alas con un al o alo del ángulo de lecha. Es o es debido a que el mallado pa a es e
ipo de alas debe ía de es a mejo adap ado y e inado. Pe o hace es os e inamien os
aumen a el iempo compu acional.
81
11. CONCLUSIONES
A la is a de los esul ados se puede conclui que el mé odo p esen ado en es e abajo supe a
ampliamen e las expec a i as iniciales. Se ha desa ollado un mé odo que consigue una
p ecisión de cálculo al a, que a oja alo es simila es a los ob enidos expe imen almen e y a los
que se pueden ob ene con o os p og amas.
El mé odo desa ollado p esen a una g an obus ez, pe mi iendo el análisis de alas de cualquie
geome ía y en ambos egímenes: es aciona io y no es aciona io.
También se han ealizado p uebas cuando el ala se some e a dis in os mo imien os en pleno
uelo. Po ejemplo las siguien es si uaciones ísicas.
- Ala a anca y pasa de es a en eposo a un ángulo de a aque 𝛼0, a mo e se con una
elocidad 𝑈∞. En un ins an e 0, el ángulo de a aque pasa a se ins an áneamen e el
doble. (Compo amien o an e un escalón)
- Ala uela a una elocidad 𝑈∞y iene además un mo imien o pe iódico en la di ección
e ical Z( ). Angulo de a aque pe manece nulo.
Pa a el caso de espues a an e un escalón en un ins an e de e minado, en el cual el ángulo de
a aque pasa de ale α=5º a un alo de α=10º, se ha ob enido la siguien e e olución de la
pendien e de la cu a de sus en ación.
FIGURA 32: RESPUESTA ANTE UN ESCALÓN
Se obse a que as la pe u bación se p oduce un descenso acusado del alo de la pendien e
de la cu a de sus en ación. Es o es lógico debido a que el ángulo de a aque pasa a se el doble
de o ma ins an ánea. La azón de que la pendien e de la cu a de sus en ación no sea
exac amen e la mi ad del alo pa a el ala en lecha, es po que el coe icien e de sus en ación
aumen a, pe o lo hace e iden emen e de o ma menos ab up a que el ángulo de a aque.
AR=4
E=1
Ψ=0
82
T as la pe u bación se puede e como el alo de la pendien e de la cu a de sus en ación
uel e a ende al alo que end ía en égimen es aciona io como e a de espe a (se puede
comp oba que el alo lími e coincide con el de la abla, ya que se ecue da que es e alo es
cons an e en égimen es aciona io según la eo ía po encial linealizada.
En cuan o a la e olución pa a el caso de mo imien o oscila o io a ángulo de a aque dado (5º)
se ha ob enido la siguien e g á ica
FIGURA 33: EVOLUCIÓN DEL COEFICIENTE DE SUSTENTACIÓN EN UN MOVIMIENTO OSCILATORIO
Donde: Ω= 𝑤𝑐
2𝑈∞
Es la ecuencia adimensional y 𝑤 la ecuencia dimensional.
Es as p uebas ya no se han podido es udia a ondo po al a de iempo y po que quedaban
ue a del obje i o p incipal del abajo. Se deja pa a es udios pos e io es la alidación de es os
esul ados.
Las condiciones mos adas an e io men e son las que se mues an en la igu a siguien e, la
cual se puede encon a en Ka z y Plo kin, 1991 [2]
FIGURA 34: EVOLUCIÓN DEL COEFICIENTE DE SUSTENTACIÓN EN UN MOVIMIENTO OSCILATORIO. (KYP) [2]
AR=4
E=1
Ψ=0
H0=0.1
Ω=0.5
α =-5º
83
12. DESARROLLOS FUTUROS
T as es e p oyec o se p oponen los siguien es desa ollos:
- Validación del mé odo con esul ados ob enidos de la li e a u a pa a el caso de
mo imien os e icales Z( ) y de cabeceo 𝛼0( ) y es udio del caso de de lexión de
supe icies hipe sus en ado as simples como laps o sla s
- Implemen ación de las ecuaciones de la elas icidad pa a esol e el p oblema
ae oelás ico pa a alas con dis in a geome ía y calcula elocidades de Flu e y de
di e gencia. Es udio del e ec o de la lecha en dichas elocidades.
- Op imización de la posición de los ale ones pa a ob ene un momen o de balance
deseado den o de los lími es de la no ma i a pa a alas de dis in a geome ía.
- Implemen ación de las ecuaciones de la elas icidad pa a esol e el p oblema
ae oelás ico pa a alas con dis in a geome ía y calcula elocidades de in e sión de
mando.
84
85
REFERENCIAS
[1]
José Manuel Go dillo A ias de Saa ed a y Guillaume Riboux Ache , "In oducción a la
Ae odinámica Po encial" Pa anin o,2012.
[2]
Joseph Ka z y Allen Plo kin, “Low-Speed Ae odynamics: F om Wing Theo y o Panel
Me hods”, McG aw-Hill, 1991.
[3]
Dowell, E. H., E. F. C awley, H. C. Cu iss, J ., D. A. Pe e s, R. H. Scanlan, and F. Sis o,
“A Mode n Cou se in Ae oelas ici y”, 3 d ed., Kluwe Academic Publishe s, 1995.
[4]
Flo es Caballe o, Míguel Ángel, “Aplicación de la Teo ía Po encial Linealizada pa a el
es udio de las ue zas eje cidas po una co ien e inciden e sob e placas ec angula es
oscilan es”, Se illa 2016
86
87
ANEXOS
1. CÓDIGO PRINCIPAL DE CÁLCULO PARA EL CASO NO ESTACIONARIO DE ALA EN
FLECHA.
%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% CÁLCULO DE COEFICIENTE DE SUSTENTACIÓN PARA ALAS EN FLECHA %%%%%%%%%%%
%% %%%%%%%%% RÉGIMEN NO ESTACIONARIO %%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% AUTOR: Fe nando Mo eno Pino %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% TUTOR : José Manuel Go dillo %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clea all; close all; clc
%% %%% INTRODUCCIÓN DE LOS PARÁMETROS
alpha0 = 5*pi/180;
ecx=[-0.238619186083197; -0.661209386466265; -0.932469514203152;0.238619186083197;
0.661209386466265; 0.932469514203152];
ecw=[0.467913934572691; 0.360761573048139; 0.171324492379170; 0.467913934572691;
0.360761573048139; 0.171324492379170];
dospi=1/2/pi;
o i=1:6
ec he ak(i)=0.5*pi*( ecx(i)+1);
end
%% GEOMETRÍA DEFINICION DE NUMERO DE PANELES Y DE INCOGNITAS
AR=5; % Ala gamien o
S=AR; % Supe icie
E=1; % Es echamien o
b=sq (S*AR); % en e gadu a
psi = 0*pi/180; % Angulo de lecha
c = 2*sq (S/AR)/(1+E); % Cue da
iala = 15; % Nume o de líneas a lo la go del eje x en el ala
Nyincog = 4*b+1; % Nume o de incógni as a lo la go del eje y
Nyincog2 = (Nyincog+1)/2;
c =c*E; % Cue da en los ex emos
i in = iala+1; % P ime a di ision de la es ela
d he ax = pi/(2*(iala-1)); % Di isión de angulo que ba e las coo denadas X
d he ay = pi/(Nyincog+1); % Di isión de ángulo que ba e las coo denadas Y
d he ay2= pi/(Nyincog2);
del ay = b/(Nyincog+1); % Sepa ación en el mallado
del ax = c/(iala-1);
Nxec = iala-2;
Nyec = Nyincog;
88
NTec = Nxec*Nyec;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% GEOMETRÍA. DEFINICION DE LA FORMA EN PLANTA Y MALLADO %%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% MALLADO DOBLE EN Y
y(1:Nyincog2+1)=(0.25*b+0.25*b*cos(((1:Nyincog2+1)-1)*d he ay2));
y(Nyincog2+1:Nyincog+2)=(-0.25*b+0.25*b*cos(((1:Nyincog2+1)-1)*d he ay2));
xa = abs(y)'* an(psi); % pun os del bo de de a aque
xs = c+(b/2* an(psi)+c -c)/(b/2)*abs(y)'; % bo de de salida
i in = iala+1; % Final de la es ela
Nxec = iala-2;
Nyec = Nyincog;
NTec = Nxec*Nyec; % El nume o o al de ecuaciones donde impone la elocidad
e ical
mx = ze os (Nyincog+2, iala+1); % Ma iz que almacena á odas las componen es x de los pun os del
mallado
my = ze os (Nyincog+2, iala+1); % Ma iz que almacena á odas las componen es y de los pun os del
mallado
UNSTEADY=1;
i UNSTEADY ==0
ies ela=1;
inT=0;
D =20*b/ies ela;
else
D =0.25;
ies ela=200; % longi ud de la es ela
inT=ies ela;
end
conT=0;
o j=1:iala
my(:,j)=y';
end
o j= 1:Nyincog+2
o i=1:iala
mx(j,i) = xa(j)+(xs(j)-xa(j))*(1-cos((i-1)*d he ax));
end
end
o j=1:Nyincog+2
o i=1:ies ela-1
i UNSTEADY ==0;
ecy(j) = my(Nyincog+3-j, iala);
89
mx (j,iala+1)=mx(j,iala)+100*b;
my(j,iala+1) = (0.5*b*cos((j-1)*d he ay));
else
mx (j,iala+i)=mx(j,iala)+i*D ;
my(j,iala+i) = my(j,iala);
end
end
end
mxS=mx(:,iala:iala+1); % P ime a anja de es ela en x
myS=my(:,iala:iala+1); % P ime a anja de es ela en y
mxE=mx(:,iala+1:end); %ma ices de la es ela en x
myE=my(:,iala+1:end); % ma ices de la es ela en y
Dxbs=mx(1,iala)-mx(1,iala-1);
%% PINTAR ALA %%
igu e (10)
o i = 1:iala % Es e bucle se enca ga de pin a las lineas en e ical de los paneles
plo ( mx(:,i), my(:, i),' ')
hold on
end
o j = 1:Nyincog+2 % Es e bucle pin a las lineas ho izon ales de los paneles
plo (mx(j,1:iala), my(j,1:iala))
hold on
end
%% PINTAR ESTELA %%
igu e(11)
o i = iala: ies ela+iala-1 % Es e bucle se enca ga de pin a las lineas en e ical de los paneles
plo ( mx(:,i), my(:, i),' ')
hold on
end
o j = 1:Nyincog+2 % Es e bucle pin a las lineas ho izon ales de los paneles
plo (mx(j,iala:ies ela+iala-1 ), my(j,iala:ies ela+iala-1))
hold on
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% CREACION DE MATRICES %%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%
sid = ((iala-2)*(Nyincog)+Nyincog); %Dimensiones de la ma iz de coe icien es en ala
%(en "sid" pun os se impone la condicion de con o no de
%impene abilidad)
sidE=((iala-2)*(Nyincog)); % Dimensión de la ma iz de coe icien es en es ela
wEIm1=ze os(sid,1); % Velocidad e ical c eada po es ela
% debido a o bellinos inyec ados en la es ela has a el ins an e -1
o j=1:Nyincog
96
ila(1) = j;
columna(1)=i;
ila(2)=j-1;
columna(2)=i;
ila(3)=j-1;
columna(3)=i-1;
ila(4)=j;
columna(4)=i-1;
end
end
o iangulo=1:2
x1=mx( ila(1), columna(1));
y1=my( ila(1), columna(1));
x2=mx( ila( iangulo+1), columna( iangulo+1));
y2=my( ila( iangulo+1), columna( iangulo+1));
x3=mx( ila( iangulo+2), columna( iangulo+2));
y3=my( ila( iangulo+2), columna( iangulo+2));
i i<=iala %No es oy en la es ela
columnaphi1= columna(1)-1+( ila(1)-2)*(iala-1);
columnaphi2= columna( iangulo+1)-1+( ila( iangulo+1)-2)*(iala-1);
columnaphi3= columna( iangulo+2)-1+( ila( iangulo+2)-2)*(iala-1);
end
D=(x2-x1)*(y3-y1)-(y2-y1)*(x3-x1);
c1=(y2-y3)/D;
d1=(x3-x2)/D;
c2=(y3-y1)/D;
d2=(x1-x3)/D;
c3 =(y1-y2)/D;
d3 =(x2-x1)/D;
s
acphi1=dospi*( I(1)+c1* I(3)+d1* I(2));
acphi2=dospi*(c2* I(3)+d2* I(2));
acphi3=dospi*(c3* I(3)+d3* I(2));
i ila(1)>1 && ila(1)<Nyincog+2 && columna(1)>1
Ma ( 0,columnaphi1)=Ma ( 0,columnaphi1)+ acphi1;
End
i ila( iangulo+1)>1 && ila( iangulo+1)<Nyincog+2 && columna( iangulo+1)>1
Ma ( 0,columnaphi2)=Ma ( 0,columnaphi2)+ acphi2;
End
i ila( iangulo+2)>1 && ila( iangulo+2)<Nyincog+2 && columna( iangulo+2)>1
Ma ( 0,columnaphi3)=Ma ( 0,columnaphi3)+ acphi3;
end
end
97
3. Mon aje_SOLAPE
%% FUNCION el mon aje en el bo de de salida y comienzo de la es ela
x0=mx(l, k); y0=my(l,k); % Coo denadas de la singula idad
0 =k-1+(l-2)*(iala-2);
Nyincog2=(Nyincog-1)/2;
ila=ze os(1,4);
columna=ze os(1,4);
i j<Nyincog2+3
i mod(i,2)==0
ipo ec angulo=1;
ila(1) = j;
columna(1)=i;
ila(2)=j-1;
columna(2)=i;
ila(3)=j-1;
columna(3)=i-1;
ila(4)=j;
columna(4)=i-1;
else
ipo ec angulo=2;
ila(4) = j;
columna(4)=i;
ila(1)=j-1;
columna(1)=i;
ila(2)=j-1;
columna(2)=i-1;
ila(3)=j;
columna(3)=i-1;
end
else
i mod(i,2)==0
ipo ec angulo=2;
ila(4) = j;
columna(4)=i;
ila(1)=j-1;
columna(1)=i;
ila(2)=j-1;
columna(2)=i-1;
ila(3)=j;
columna(3)=i-1;
else
ipo ec angulo=1;
ila(1) = j;
columna(1)=i;
98
ila(2)=j-1;
columna(2)=i;
ila(3)=j-1;
columna(3)=i-1;
ila(4)=j;
columna(4)=i-1;
end
end
o iangulo=1:2
x1=mxS( ila(1), columna(1));
y1=myS( ila(1), columna(1));
x2=mxS( ila( iangulo+1), columna( iangulo+1));
y2=myS( ila( iangulo+1), columna( iangulo+1));
x3=mxS( ila( iangulo+2), columna( iangulo+2));
y3=myS( ila( iangulo+2), columna( iangulo+2));
i i<=2 %No es oy en la es ela
columnaphi1= columna(1)+( ila(1)-2)*(2);
columnaphi2= columna( iangulo+1)+( ila( iangulo+1)-2)*(2);
columnaphi3= columna( iangulo+2)+( ila( iangulo+2)-2)*(2);
end
D=(x2-x1)*(y3-y1)-(y2-y1)*(x3-x1);
c1=(y2-y3)/D;
d1=(x3-x2)/D;
c2=(y3-y1)/D;
d2=(x1-x3)/D;
c3 =(y1-y2)/D;
d3 =(x2-x1)/D;
s
acphi1=dospi*( I(1)+c1* I(3)+d1* I(2));
acphi2=dospi*(c2* I(3)+d2* I(2));
acphi3=dospi*(c3* I(3)+d3* I(2));
i ila(1)>1 && ila(1)<Nyincog+2 && columna(1)>1
Ma SOLAPE( 0,columnaphi1)=Ma SOLAPE( 0,columnaphi1)+ acphi1;
End
i ila( iangulo+1)>1 && ila( iangulo+1)<Nyincog+2 && columna( iangulo+1)>0
Ma SOLAPE( 0,columnaphi2)=Ma SOLAPE( 0,columnaphi2)+ acphi2;
End
i ila( iangulo+2)>1 && ila( iangulo+2)<Nyincog+2 && columna( iangulo+2)>0
Ma SOLAPE( 0,columnaphi3)=Ma SOLAPE( 0,columnaphi3)+ acphi3;
end
end
99
4. Mon aje_Es ela_2
%% FUNCION PARA COLOCAR EN LA MATRIZ DE LA ESTELA
x0=mx(l, k); y0=my(l,k); % Coo denadas de la singula idad
0 =k-1+(l-2)*(iala-2);
Nyincog2=(Nyincog-1)/2;
ila=ze os(1,4);
columna=ze os(1,4);
i j<Nyincog2+3
i mod(i,2)==0
ipo ec angulo=1;
ila(1) = j;
columna(1)=i;
ila(2)=j-1;
columna(2)=i;
ila(3)=j-1;
columna(3)=i-1;
ila(4)=j;
columna(4)=i-1;
else
ipo ec angulo=2;
ila(4) = j;
columna(4)=i;
ila(1)=j-1;
columna(1)=i;
ila(2)=j-1;
columna(2)=i-1;
ila(3)=j;
columna(3)=i-1;
end
else
i mod(i,2)==0
ipo ec angulo=2;
ila(4) = j;
columna(4)=i;
ila(1)=j-1;
columna(1)=i;
ila(2)=j-1;
columna(2)=i-1;
ila(3)=j;
columna(3)=i-1;
else
ipo ec angulo=1;
ila(1) = j;
columna(1)=i;
100
ila(2)=j-1;
columna(2)=i;
ila(3)=j-1;
columna(3)=i-1;
ila(4)=j;
columna(4)=i-1;
end
end
o iangulo=1:2
x1=mxE( ila(1), columna(1));
y1=myE( ila(1), columna(1));
x2=mxE( ila( iangulo+1), columna( iangulo+1));
y2=myE( ila( iangulo+1), columna( iangulo+1));
x3=mxE( ila( iangulo+2), columna( iangulo+2));
y3=myE( ila( iangulo+2), columna( iangulo+2));
i i<=ies ela-1 %No es oy en la es ela
columnaphi1= columna(1)+( ila(1)-2)*(ies ela-1);
columnaphi2= columna( iangulo+1)+( ila( iangulo+1)-2)*(ies ela-1);
columnaphi3= columna( iangulo+2)+( ila( iangulo+2)-2)*(ies ela-1);
end
D=(x2-x1)*(y3-y1)-(y2-y1)*(x3-x1);
c1=(y2-y3)/D;
d1=(x3-x2)/D;
c2=(y3-y1)/D;
d2=(x1-x3)/D;
c3 =(y1-y2)/D;
d3 =(x2-x1)/D;
s
acphi1=dospi*( I(1)+c1* I(3)+d1* I(2));
acphi2=dospi*(c2* I(3)+d2* I(2));
acphi3=dospi*(c3* I(3)+d3* I(2));
i ila(1)>1 && ila(1)<Nyincog+2 && columna(1)>1
Ma E( 0,columnaphi1)=Ma E( 0,columnaphi1)+ acphi1;
End
i ila( iangulo+1)>1 && ila( iangulo+1)<Nyincog+2 && columna( iangulo+1)>0
Ma E( 0,columnaphi2)=Ma E( 0,columnaphi2)+ acphi2;
end
i ila( iangulo+2)>1 && ila( iangulo+2)<Nyincog+2 && columna( iangulo+2)>0
Ma E( 0,columnaphi3)=Ma E( 0,columnaphi3)+ acphi3;
end
end
101
5. Mon aje_0
%% FUNCION QUE COLOCA LOS VALORES DE LA MATRIZ CUANDO ANALIZAMOS LA SINGULARIDAD
x0=mx(l,k);
y0=my(l,k);
x1=x0;
y1=y0;
ila(1)=l;
columna(1)=k;
ila(2)=l;
columna(2)=k+1;
ila(3)=l-1;
columna(3)=k+1;
ila(4)=l-1;
columna(4)=k;
ila(5)=l-1;
columna(5)=k-1;
ila(6)=l;
columna(6)=k-1;
ila(7)=l+1;
columna(7)=k-1;
ila(8)=l+1;
columna(8)=k;
ila(9)=l+1;
columna(9)=k+1;
ila(10)=l;
columna(10)=k+1;
ila0=k-1+(l-2)*(iala-2);
columna0=k-1+(l-2)*(iala-1);
o iangulo=1:8
x2=mx( ila( iangulo+1), columna( iangulo+1));
y2=my( ila( iangulo+1), columna( iangulo+1));
x3=mx( ila( iangulo+2), columna( iangulo+2));
y3=my( ila( iangulo+2), columna( iangulo+2));
columnaphi2= columna( iangulo+1)-1+( ila( iangulo+1)-2)*(iala-1);
columnaphi3= columna( iangulo+2)-1+( ila( iangulo+2)-2)*(iala-1);
D=(x2-x1)*(y3-y1)-(y2-y1)*(x3-x1);
c1=(y2-y3)/D;
d1=(x3-x2)/D;
c2=(y3-y1)/D;
d2=(x1-x3)/D;
c3=(y1-y2)/D;
d3=(x2-x1)/D;
s_0
acphi1=dospi*( I(1));
Ma ( ila0,columna0)=Ma ( ila0,columna0)+ acphi1;
end
102
6. FS
%% PROGRAMA QUE REALIZA LA INTEGRACIÓN EN CASO DE ESTAR EN UN PUNTO QUE NO ES LA
%%SINGULARIDAD %%%
1 0x=x1-x0;
1 0y=y1-y0;
2 1x=x2-x1;
2 1y=y2-y1;
3 1x=x3-x1;
3 1y=y3-y1;
3 2x=x3-x2;
3 2y=y3-y2;
L12=sq ( 2 1x* 2 1x+ 2 1y* 2 1y);
L10=sq ( 1 0x* 1 0x+ 1 0y* 1 0y);
he a2=angulo( 2 1x, 2 1y);
he a3=angulo( 3 1x, 3 1y);
i (abs( he a3)<1e-12 && he a3< he a2)
he a3=2*pi;
end
he a4=angulo( 3 2x, 3 2y);
i he a2<0|| he a3<0|| he a4<0
disp('pa a')
pause
end
I1k=0;
I2k=0;
I3k=0;
o p=1:6
he ak=0.5*( he a3- he a2)*( ecx(p)+1)+ he a2;
k=sin( he a4- he a2)*L12/(sin( he a4)*cos( he ak)-cos( he a4)*sin( he ak));
x0k=cos( he ak)* 1 0x+sin( he ak)* 1 0y;
xk= k+ x0k;
a2k=L10*L10- x0k* x0k;
S0k=1/sq ( xk* xk+ a2k)-1/L10;
i (abs( a2k)>1e-8)
S3k=(1/ a2k)*( xk/sq ( a2k+ xk* xk)- x0k/sq ( a2k+ x0k* x0k));
S4k=log((sq ( a2k+ xk* xk)+ xk)/(sq ( a2k+ x0k* x0k)+ x0k));
S1k=-S0k- x0k*S3k;
S2k=S4k-2* x0k*(-S0k- x0k*S3k)-L10*L10*S3k;
else
S3k=-0.5*(1/(( k+L10)*( k+L10))-1/(L10*L10));
S4k=log((L10+ k)/L10);
S1k=-S0k- x0k*S3k;
S2k=S4k-2* x0k*(-S0k- x0k*S3k)-L10*L10*S3k;
103
end
I1k= I1k+S1k* ecw(p);
I2k= I2k+S2k*sin( he ak)* ecw(p);
I3k= I3k+S2k*cos( he ak)* ecw(p);
end
I(1)=0.5*( he a3- he a2)* I1k;
I(2)=0.5*( he a3- he a2)* I2k;
I(3)=0.5*( he a3- he a2)* I3k;
104
7. FS_0
%% INTEGRACION EN CASO DE SINGULARIDAD
2 1x= x2-x1;
2 1y= y2-y1;
3 1x= x3-x1;
3 1y=y3-y1;
3 2x=x3-x2;
3 2y=y3-y2;
L12 = sq ( 2 1x* 2 1x+ 2 1y* 2 1y);
he a2=angulo( 2 1x, 2 1y);
he a3=angulo( 3 1x, 3 1y);
i (abs( he a3)<1e-12 && he a3< he a2)
he a3=2*pi;
end
he a4=angulo( 3 2x, 3 2y);
I1k=0;
o p=1:6
he ak=0.5*( he a3- he a2)*( ecx(p)+1)+ he a2;
k=sin( he a4- he a2)*L12/(sin( he a4)*cos( he ak)-cos( he a4)*sin( he ak));
I1k= I1k+(-1/ k)* ecw(p);
end
I(1)=0.5*( he a3- he a2)* I1k;
105
8. Ángulo
%% FUNCION QUE COGE EL ÁNGULO CORRECTO PARA LA INTEGRACION %%
unc ion [ he a]=angulo( x, y)
he a0=a an(abs( y)/abs( x));
i (abs( x)>1e-15)
he a0=a an(abs( y)/abs( x));
else
i ( y>0)
he a=pi/2.;
else
he a=3*pi/2.;
end
end
i (abs( y)<1e-15)
i ( x>0)
he a=0.;
else
he a=pi;
end
else
i ( y>1e-15 && x>1e-15)
he a= he a0;
end
i ( y>1e-15 && x<-1e-15)
he a=pi- he a0;
end
i ( y<1e-15 && x<-1e-15)
he a=pi+ he a0;
end
i ( y<1e-15 && x>1e-15)
he a=2*pi- he a0;
end
end
end
112
acphi2=dospi*(c2* I(3)+d2* I(2));
acphi3=dospi*(c3* I(3)+d3* I(2));
i ila(1)>1 && ila(1)<Nyincog+2 && columna(1)>1
Ma ( 0,columnaphi1)=Ma ( 0,columnaphi1)+ acphi1;
End
i ila( iangulo+1)>1 && ila( iangulo+1)<Nyincog+2 && columna( iangulo+1)>1
Ma ( 0,columnaphi2)=Ma ( 0,columnaphi2)+ acphi2;
End
i ila( iangulo+2)>1 && ila( iangulo+2)<Nyincog+2 && columna( iangulo+2)>1
Ma ( 0,columnaphi3)=Ma ( 0,columnaphi3)+ acphi3;
end
end