MM

Semestrální práce z předmětu

KMA/MM

Jméno: Plánička Stanislav Mail: [email protected] Os. Číslo: A08N0144P

Obsah OBSAH

2

FYZIKÁLNÍ MINIMUM

3

KINEMATIKA Idealizace hmotným bodem Zlatá rovnice kinematiky DYNAMIKA 2. Newtonův pohybový zákon AERODYNAMIKA Základ mechaniky tekutin Aerodynamické síly působící na křídlo Mezinárodní standardní atmosféra MSA ŘEŠENÍ – APLIKACE ZNALOSTÍ NA LETOUN Hustota vzduchu, vliv gravitace Vnější síly působící na letoun Sestavení pohybové rovnice a výpočet kinematických veličin ŘEŠENÍ A SIMULACE LETU LETADLA Vstupní parametry Určení aerodynamických součinitelů Simulace letu letadla MODIFIKOVANÝ MODEL Přehled výkonů letadla

3 3 3 4 4 4 4 5 8 8 9 10 11 12 12 13 13 16 18

ZÁVĚR

19

POUŽITÁ LITERATURA

20

2

Fyzikální minimum Kinematika Idealizace hmotným bodem S pojmem „hmotný bod“ se ve fyzice a mechanice setkáváme poměrně často. Jedná se o idealizovaný objekt, který má pouze hmotnost soustředěnou do daného bodu. Jeho vnější rozměry jsou nulové. Tato idealizace má svoje opodstatnění. Hmotným bodem nahrazujeme těleso, jehož rozměry a tvar jsou pro řešení sledovaného pohybu nepodstatné. Velmi užitečná je zvláště v případech, kdy těleso koná „posuvný pohyb“. Polohu libovolného bodu daného tělesa lze zapsat jako vektorový součet vektoru polohy referenčního bodu (je užitečné aby tento bod byl zároveň těžištěm) a polohy tohoto bodu vůči referenčnímu r r r R = r + r *. (1)

Obr.1: Polohu libovolného bodu tělesa lze jednoznačně určit prostřednictvím polohy referenčního bodu a vzájemné polohy obou.

r

Koná –li těleso posuvný pohyb, je relativní vektor r * v čase konstantní. Dá se snadno ukázat, že v takovém případě opisují všechny body tělesa stejné, vzájemně posunuté trajektorie. Jejich rychlosti a zrychlení (kinematické poměry) jsou totožné – lze se o tom přesvědčit prostým derivováním rovnosti (1). V případě posuvného pohybu tedy neztrácíme (z pohledu kinematiky) redukcí tělesa na hmotný bod žádnou podstatnou informaci.

Zlatá rovnice kinematiky Vztahy mezi kinematickými veličinami si všichni jistě pamatujeme z hodin fyziky. r v Připomeňme, že platí: vektor rychlosti v je časovou derivací polohového vektoru r r r d r r& v= =r, (2) dt r vektor zrychlení a je časovou derivací vektoru rychlosti r r d v r& a= =v, (3) dt tyto rovnosti mechanici nazývají „zlatá rovnice kinematiky“ a zapisují: r r r d v d2 r a= = . (4) dt dt2 V kartézském souřadnicovém systému tyto vztahy s výhodou rozepisujeme po složkách.

3

Dynamika 2. Newtonův pohybový zákon „Zákon síly“ je snad nejvýznamnější přírodní princip používaný v mechanice. Jeho diferenciální tvar formulujeme matematicky: r dp r =F, (5) dt jedná se o lokální zápis zákona bilance hybnosti. Pokud je pravá strana nulová (tzv. izolované soustavy) obdržíme matematický zápis zákona zachování hybnosti – jednoho z elementárních přírodních principů. Zákon síly je pro nás důležitý, neboť nám umožňuje přechod od dynamiky ke kinematice, r jak nyní ukážeme. Hybnost p je definována vztahem r r p = mv ; (6) m značí hmotnost, kterou v klasické mechanice tuhého tělesa považujeme za konstantní. Dosazením (6) do (5) obdržíme zákon síly ve známém tvaru r r F = ma . (7) r r Tato jednoduchá formule dává do souvislosti F a zrychlení a jednoduchým vynásobením r konstantou m . F je výsledná síla, která na těleso působí, je to vektorový součet všech r působících sil Fi (jejich počet označíme n ) n r r (8) F = ∑ Fi . i =1

V případech, kdy je působící silová soustava komplikovanější, s výhodou použijeme jinou r formulaci téhož tzv. d’Alambertův princip, který zavádí „setrvačnou sílu“ S : r r S = −m a (9) D’Alambertův zákon tvrdí, že součet působících sil a síly setrvačné je nulový: r r r F+S =0 (10)

Aerodynamika Základ mechaniky tekutin Abychom mohli porozumět vzniku aerodynamických sil působících na letadlo od obtékajícího proudu, je třeba mít správnou fyzikální představu o proudovém poli a rozložení tlaků v jeho blízkosti. K tomu nám poslouží dva jednoduché vztahy, které jsou, za jistých zjednodušujících podmínek, odvozeny z fundamentálních rovnic mechaniky tekutin tzv. Naviérových – Stokesových rovnic. Prvním vztahem, který představuje lokální formulaci zákona zachování hmoty, je „rovnice kontinuity“ ∂ρ r + div(ρ v ) = 0 (11) ∂t Pro zvolenou proudovou trubici proměnného průřezu S ( x) , za předpokladu stacionárního proudění a nestlačitelné tekutiny ji můžeme upravit do tvaru S v = konst. (12) Zápis rovnice kontinuity v tomto tvaru ukazuje důležitý závěr: v místě zúžení proudové trubice je větší rychlost a naopak v místech širších rychlost klesá. Poznámka: Pokud jsme pozorní zajisté jsme se pozastavili nad nestlačitelností tekutiny. Tekutinou rozumíme v aerodynamice vzduch, který je poměrně snadno stlačitelný. Jeho stlačitelnost ovšem závisí na stavu určeném

4

především „Machovým číslem“ rychlostí zvuku c

M , definovaným jako poměr mezi velikostí rychlosti proudu v a místní

M =

v . c

(13)

Pokud je hodnota Machova čísla malá M < 0,1 −& 0,3 lze považovat vzduch s dostatečnou přesností za nestlačitelné médium. Tato podmínka je splněna u většiny běžných menších a sportovních letadel.

Tím druhým významným vztahem je „Bernoulliho rovnice“ pro proudovou trubici, která vychází z již zmíněného zákona síly (potažmo pohybových Navierových-Stokesových rovnic) s přijetím některých zjednodušujících předpokladů (proudění je stacionární, tekutina nevazká atd.). My vyjdeme ze vztahu 1 p + ρ v 2 = konst. , (14) 2 ve kterém první člen nazýváme „statický tlak“, druhý „dynamický tlak“. Na základě analýzy vztahu můžeme tedy tvrdit, že v místě s vyšší rychlostí proudu zaznamenáme pokles tlaku, naopak v místě zúžení statický tlak poroste.

Aerodynamické síly působící na křídlo Každé těleso, pohybující se v určitém prostředí, je vystaveno působení sil vznikajících při jeho pohybu. Přitom nezáleží na tom, zda je nehybný předmět obtékán prostředím, nebo zda se pohybuje těleso v klidném prostředí (rychlostí opačného směru a stejné velikosti). Při obtékání těles symetrických ve směru proudu působí výsledná síla vzduchu, kterou nazýváme r „výsledná aerodynamická síla“ R , ve směru proudu. Při obtékání těles nesymetrických (popř. nesymetrických vůči smyslu proudu) nepůsobí již výsledná aerodynamická síla ve r směru proudu, ale je od něj odchýlena. V tomto případě rozkládáme sílu R na dvě složky: r složku působící ve směru proudu nazýváme „aerodynamický odpor“ a značíme D , a složku r kolmou vůči směru proudu „aerodynamický vztlak“ L ; platí tedy r r r L+D= R (15) Aerodynamické síly působící na letadlo zkoumá aerodynamika letadla - obor, který se vyčlenil z mechaniky tekutin. Jejím úkolem je řešit proudové pole v okolí letadla, stanovit síly působící na jeho povrch (úlohy přímé), či na základě požadovaných sil, nebo proudového pole stanovit optimální tvar (nepřímá úloha). V běžných případech bývá obtékajícím prostředím vzduch zemské atmosféry. Z výsledků experimentů a numerických simulací víme, že při rovinném laminárním obtékání profilu křídla vznikne charakteristické proudové pole.

Obr.2: Rozložení proudnic při rovinném laminárním obtékání křídla letadla. Povšimněte si, že na horní straně křídla je větší hustota proudnic, proudové trubice jsou r zúženy (oproti nerozrušenému nabíhajícímu proudu o rychlosti v∞ ). Na spodní straně křídla jsou proudové trubice mírně rozšířené. Zúžení trubice znamená – viz rovnice kontinuity – vyšší rychlost a tedy – podle Bernoulliho rovnice – nižší tlak. Analogickou úvahou usuzujeme

5

na přetlak působící ne spodní stranu křídla. Tento závěr je zcela zásadní: Při ofukování profilu proudem vzduchu vzniká na jeho horní straně oblast sníženého statického tlaku (podtlaku), na spodní straně oblast zvýšeného statického tlaku (přetlaku). Letící křídlo je tedy shora vzduchem přisáváno a zespodu tlačeno. r Vzduch obtékající letadlo rychlosti v∞ (charakterizovaný hustotou ρ ∞ , tlakem p∞ , termodynamickou teplotou T∞ a za pohybu dynamickou viskozitou η ) vytvoří při obtékání křídla na jeho povrchu proměnné rozložení statického tlaku p a smykového napětí τ . (Povšimněme si, že v aplikované aerodynamice zavedená terminologie odpovídá nehybnému r křídlu a proudu vzduchu, který jej obtéká rychlostí v∞ ).

Obr.3: Rozložení tlaků působících na plochu křídla. Z rozložení tlaků na profilu vidíme, že přibližně dvě třetiny celkového rozdílu generuje podtlak nad křídlem. Sání na horní straně profilu má významný podíl na vztlaku. r r Na nekonečně malou plošku letadla o velikosti dS působí elementární síly dFp a dF f . První zveme tlakovou a druhá zvaná smyková (třecí). Integrací obou těchto složek po celém r r povrchu dostaneme výslednou aerodynamickou sílu R . Vektor síly R tedy obvykle r r rozkládáme do složek D ve směru nabíhajícího proudu a do složky k ní kolmé L . (Výsledná síla se alternativně rozkládá, především pro použití v pružnostářských a pevnostních výpočtech, ve směru letadla, popř. tětivy křídla a směru kolmého.) Obě tyto síly vypočítáváme empirickými vzorci: 1 1 2 2 L = c L ρ S v∞ a D = c D ρ S v∞ , (16, 16a) 2 2 kde c L a c D jsou bezrozměrné součinitele vztlaku a odporu. S je nosná plocha. Při zvyšování rychlosti dochází k tomu, že proudění vzduchu v okolí letadla ztrácí svůj hladký, laminární charakter a vyvíjí se proudění turbulentní. Tyto nepravidelnosti v proudovém poli, spojené s vířením a úplavovýmy jevy znamenají další růst odporových sil. Roste – li rychlost letadla dále, nelze již přijmout platnost modelu s nestlačitelnou tekutinou. Matematický model popisující proudové pole se ovšem zahrnutím stlačitelnosti nezanedbatelně zesložiťuje. V subsonické oblasti ( M < 1 ) se charakter proudu zásadně (kvalitativně) neodlišuje od proudění s nestlačitelnou tekutinou. Ovšem při přiblížení k této hranici dochází k bouřlivým změnám. Charakter proudění se významně mění, vznikají rázové vlny, které zvyšují odpor a negativně působí na vztlak. (Tyto děje jsou důsledkem právě stlačitelnosti!) Krom toho již v subsonických oblastech dochází k tzv. odtržení proudu, které 6

opět snižuje aerodynamickou účinnost křídla. Tento jev je značně závislý na tvaru křídla a režimu letu letadla. Obecný matematický popis těchto jevů je tedy značně problematický. Rovnice pro výpočet odporu a vztlaku by patrně byly pouze přiblížením silně závislým na případě. Z výše uvedeného můžeme usuzovat na jejich nemonotónnost, nehladkost a dokonce nespojitost. Bohužel žádná dostupná literatura se nevěnuje ani nejmenšímu rozšíření Newtonových vztahů (16) do transsonické oblasti ( M > 1 ), všechny publikace předpokládají jeho platnost v celém rozsahu rychlostí. I nám tedy nezbude, než se spokojit s tímto empirickým přiblížením.

Obr.4: Rozložení statického tlaku v okolí obtékaného křídla – porovnání rozdílu modelů s nestlačitelnou (vlevo) a stlačitelnou tekutinou (vpravo). Můžeme si všimnout podobného charakteru tlakového pole obou případů. Numerická simulace z CFD výpočtového systému Fluent. (Pozn.: simulace „vlevo“ je při rychlosti nabíhajícího proudu v = 100 m s-1, simulace „vpravo“ je při rychlosti nabíhajícího proudu v =& 90 m s-1)

Obr.5: Rychlostní pole v okolí obtékaného křídla při rychlosti nabíhajícího proudu v = 100 m s-1 (model s nestlačitelnou tekutinou). Numerická simulace z CFD výpočtového systému Fluent.

Obr.6: Tlakové pole při mírně nadzvukové (transsonické) rychlosti. Všimněte si tlakové vlny utvářející se před profilem a na odtokové hraně profilu. Jedná se samozřejmě o model respektující stlačitelnost vzduchu. Numerická simulace z CFD výpočtového systému Fluent.

7

Mezinárodní standardní atmosféra MSA Již krátce po letu Flyeru vyvstala potřeba jakési standardní atmosféry, pomocí které by se daly vypočítávat a porovnávat výkony letadel. Tento „etalon“ prošel pochopitelně vývojem a nyní je MSA srovnávací atmosféra dána normou ISO z roku 1975. Mezi informace standardní atmosféry je třeba zařadit především geofyzikální údaje o zeměkouli. Není to přesná koule, ale zploštělé rotační těleso zvané zemský elipsoid, nebo též geoid, jež je přibližně symetrický vůči rovině zemského rovníku. Má tvar kružnice o poloměru 6378,16 km. Jeho poloosa symetrie, jež je současně osou zemské rotace, má délku 6356,77 km. Za definiční místo normálního stavu atmosférického vzduchu bylo zvoleno místo na mořské hladině o zeměpisné šířce ϕ 0 = 45°32′33′′ , kde je střední poloměr Země r0 = 6356,766 km a tíhové zrychlení g 0 = 9,80665 m s −2 . Vertikální vzdálenost od tohoto místa se nazývá „geometrická nadmořská výška“ h . Střední hodnoty stavových veličin vzduchu v tomto místě jsou: termodynamická teplota T0 = 288,15 K , tlak p 0 = 101325 Pa , hustota ρ 0 = 1,2255 kg m −3 . Závislost gravitačního zrychlení na geometrické výšce stanovíme z Newtonova gravitačního zákona . Podle něj je gravitační síla mezi dvěma hmotnými body nepřímo úměrná kvadrátu jejich vzdáleností. Východiskem ke stanovení závislosti termodynamických stavových veličin atmosférického vzduchu na nadmořské výšce je sestavení diagramu středních hodnot teploty v závislosti na nadmořské výšce získaného z dostatečného počtu výškových měření. Kontrolní měření tlaku ukázala, že vzduch lze ve výškovém rozsahu praktického letectví pokládat za ideální plyn. Je tedy možno použít standardní stavovou rovnici p = rT ; (17)

ρ

r je měrná plynová konstanta s hodnotou pro vzduch r = 287,05287 J kg-1 K-1. Stejně tak

v těchto výškách (cca h < 25 km) lze hustotu vzduchu ve výšce h * (což je geometrická nadmořská výška v kilometrech) vypočíst pomocí polynomu třetího stupně 2 3 ρ (h ) = m + n h* + o h* + p h* ; (18)

ρ0

m = 1 , n = −0,0965 km , o = 0,00325 km-2, p = −3,7 ⋅ 10 −5 km-3. Pozor, do vztahu (18) je -1

potřeba dosadit h * v kilometrech!

Řešení – Aplikace znalostí na letoun Nyní již víme vše potřebné. Využijeme našich znalostí a pokusíme se vytvořit model popisující aerodynamiku a mechaniku letu letadla. V úvodní kapitole jsme postupovali od obecných poznatků k detailům. V této se ale nejprve zaměříme na detaily, bez nichž by nebylo možné model sestavit a od nich pak přejdeme k určení působících sil, sestavíme pohybovou rovnici a zformulujeme závislosti kinematických veličin. Pohyb budeme uvažovat pouze v rovině, tuto rovinu nazveme „výpočtová oblast“ Ω . Polohu letadla bude jednoznačně určovat kartézský souřadnicový systém 0xy. Letadlo tedy nebude zatáčet, zároveň zanedbáme vliv křivosti zemského povrchu. Dále budeme uvažovat, že tíhové pole je homogenní (ekvipotenciály jsou rovnoběžné s osou x) a termodynamické veličiny určující stav atmosféry se mění pouze ve vertikálním smyslu (pro parciální derivace všech veličin ve směru x platí ∂A / ∂x = 0 ) podle vztahů MSA.

8

Hustota vzduchu, vliv gravitace Pro sestavení pohybové rovnice musíme znát vnější síly působící na letoun. Pro správné vyčíslení aerodynamických sil budeme potřebovat pomoc MSA, neboť určení jejich hodnot dle (16) vyžaduje znalost „místní“ hustoty nabíhajícího proudu, která je silně závislá na geometrické výšce (letové hladině). Velikost hustoty vypočteme ze vztahu (18).

Obr.7: Polynomická závislost hustoty na nadmořské výšce podle MSA. Že vliv gravitace klesá s výškou nás asi nepřekvapí. Všichni víme, že kosmonauti ve velkých výškách zažívají tzv. „stav bez tíže“ (tento pojem je značně zavádějící). Námi uvažované výšky jsou však mnohem menší a tedy i pokles velikosti gravitačního zrychlení bude menší. Pro jednoduchost ztotožníme vektory gravitačního a tíhového zrychlení, jeho vektor má v každém místě oblasti Ω směr opačný s osou y zavedeného souřadnicového r r r r systému ( ∀ X ∈ Ω, ∀ t ∈ 0, T : g ( X , t ) ≡ konst. = (0,− g ), g • e y = − g = − g ). Velikost vektoru g v závislosti na výšce h určíme z Newtonova gravitačního zákona, v souladu s předpoklady přijatými v MSA viz výše. 2 r0 1 gh = g0 =κ M (19) 2 2 r0 + hgeom r0 + hgeom

(

)

(

)

Obr. 8: Závislost tíhového zrychlení na nadmořské výšce podle Newtonova gravitačního zákona (popř. MSA). Z obr. 8 nebo dosazením do rovnice (19) zjistíme, že velikost tíhy klesá ve výškovém rozsahu praktického letectví jen nepatrně. Ve výšce h = 25 km (velká výška dosahovaná specializovanými „výškovými“ vojenskými letouny) má tíhové zrychlení velikost

9

g 25 km = 9,7300 m s-2, což je rozdíl oproti referenční hodnotě g 0 menší než 1%. Nebudeme tedy nadále uvažovat závislost gravitačních sil na výšce. Vektor tíhového zrychlení v celé výpočtové oblasti Ω budeme považovat za konstantní, jeho velikost udáme na g = 9,78 m s-2, čímž vnášíme chybu v řádech X ‰, která nezhoršuje přesnost výpočtů.

Vnější síly působící na letoun Na letoun působí 3 druhy vnějších sil: aerodynamické síly (vztlak a odpor), tíhová síla a hnací síla motoru. Protože předpokládáme, že letoun koná rovinný, posuvný pohyb, využijeme rozklad těchto sil do směrů x a y.

Obr. 9: Síly působící na letoun. Zeleně jsou vyznačeny aerodynamické síly (L – vztlak, lift; D – odpor, drag), červeně tahová síla motoru (M), modře tíhová síla (G). Tlusté šedivé šipky představují setrvačné síly zavedené v opačném smyslu černých souřadnicových os x a y. Výslednice vnějších sil ve směru x je (v libovolném časovém okamžiku letu a libovolném bodě výpočtové oblasti) výsledkem součtu sil působících ve směru x: Fx = M − D , (20) analogicky výslednice sil ve směru y: Fy = −G + L , (20a) rovnice (20) a (21) již nejsou vektorové, ale skalární. „Tahová síla motoru“ M budeme v našem případě dosahovat pouze tří diskrétních hodnot M ∈ {0, M 100 , M F } . U reálných běžných letadel je tah motoru funkcí nastavení jeho výkonu: vypnutý motor pochopitelně neprodukuje tah; ve standardním pracovním režimu je tah spojitou, kladnou, shora i zdola omezenou, rostoucí funkcí výkonu; přídavné spalování tzv. forsáž (pouze některá vojenská letadla) umožňuje další navýšení nad „100% tah“. Pro jednoduchost budeme uvažovat, že se tah s narůstající výškou (při jeho neměnném pracovním režimu) nemění - v reálu klesá s rostoucí výškou letu (u letadel „poháněných“ vrtulí dosti významně). „Tíhová síla“ se (pro ∀t ∈ 0, T ) vypočte vynásobením okamžité hmotnosti ( m(t ) ) tíhovým zrychlením

10

G = mg . (21) Při letu bude hmotnost časem klesat, protože její část letoun ztrácí spalováním paliva. Rychlost tohoto úbytku (v jednotkách kg s-1) závisí na pracovním režimu motoru a lze ji vyjádřit jako časovou derivaci celkové hmotnosti dm ∆m = m& = . (22) dt V našem případě tedy bude rychlost změny hmotnosti nabývat pouze tří nekladných hodnot: ∆m ∈ {0, ∆m100 , ∆mF }. Hmotnost v čase t ∈ 0, T určíme jako počáteční hmotnost zmenšenou o celkový úbytek hmotnosti vlivem spalování v časovém intervalu τ ∈ 0, t : t

m = m0 + ∫ ∆m dτ .

(23)

0

K čemuž dojdeme též integrací rovnice (22). Poznámka: Při skládání sil v rovině je nutné tyto dvě tzv. složkové podmínky doplnit podmínkou momentovou. Kdyby přesně platilo schéma na obr. 9, byla by tato podmínka triviální. Je ale dobré si uvědomit některé podstatné detaily, které tento model opomíjí. Aerodynamické síly v našem modelu jsou síly již vzniklé složení dílčích aerodynamických sil (křídla, vodorovná ocasní plocha a trup).

L

Obr. 10: Skutečné rozložení sily působících na letadlo ve směru y. V ideálním případě je nositelka výsledného vztlaku totožná s nositelkou tíhové síly. Ale u opravdových letadel není obecně poloha těžiště prázdného letounu totožná s polohou těžiště paliva, resp. natankovaného letounu (o užitečné zátěži ani nemluvě). Během letu tedy dochází k posunu polohy těžiště. Případná neshoda poloh nositelek generuje nenulový moment. Obdobně můžeme přemýšlet o silách ve směru osy x. Je tedy nasnadě, že ve skutečnosti obecně nebude momentová podmínka triviální. Tento nesoulad se vyrovnává dalšími aerodynamickými silami na řídících plochách (klapky, křidélka, VOP, aerodynamické brzdy aj.), což je spojeno se změnou vztlaku (pozitivní i negativní) a růstem odporu (vždy). Pilot s těmito jevy musí počítat. U moderních letounů s digitálním systémem řízením po drátě (fly-by-wire) řeší tyto problémy řídící jednotka automaticky. Reálný model silového působení lze tedy převést na náš (na obr. 9), za cenu toho, že bude třeba dále opravovat hodnoty vztlaku a odporu. Ale jelikož výsledný moment nebývá velký a my budeme předpokládat vynikající práci konstruktérů letadel, nebudeme do modelu zahrnovat efekty plynoucí z momentové podmínky rovnováhy.

Sestavení pohybové rovnice a výpočet kinematických veličin Protože náš model letu letadla předpokládá vodorovný let v celém jeho trváním (což je posuvný pohyb), použijeme zjednodušení hmotným bodem. 11

Budeme uvažovat jednoduché počáteční podmínky. V čase t = 0 s platí: • poloha letounu (který je nahrazen hmotným bodem, jehož poloha je totožná s jednotným působištěm všech sil na letoun) je totožná s počátkem souřadnicového systému • hodnoty všech kinematických veličin jsou nulové Po zavedení setrvačných sil (viz obr. 9, (9)) využijeme k sestavení pohybové rovnice d’Alambertův princip (10). Skalární pohybová rovnice ve směru x: M − D − max = 0 (24) skalární pohybová rovnice ve směru y: − G + L − may = 0 (24a) Z těchto rovnic vyjádříme hledaná zrychlení, která se stanou vstupními veličinami do kinematických závislostí. M (t ) − D(v(t ), ρ (h(t ))) − G(m(t )) + L(v(t ), ρ (h(t ))) a x (t ) = , a y (t ) = (25, 25a) m(t ) m(t ) r r Při hledání hodnot kinematických veličin (rychlost v , poloha r ) dosadíme zrychlení do zlaté rovnice kinematiky, kterou rozepíšeme po složkách. Rychlosti v jednotlivých směrech v libovolném časovém okamžiku ∀t ∈ 0, T t

t

v x = ∫ a x (τ ) d τ , v y = ∫ a y (τ ) d τ 0

(26, 26a)

0

r příslušné polohy (složky polohového vektoru r ) vypočteme v ∀t ∈ 0, T : t

t

x = ∫ v x (τ ) d τ , y = ∫ v y (τ ) d τ . 0

(27, 27a)

0

Řešení soustavy rovnic (25), (26) a (27) není kvůli jejich složité „implicitní“ formulaci jednoduché, nicméně by nám poskytlo kompletní informace o okamžité poloze, rychlosti a zrychlení letounu. Proto se níže pokusíme alespoň o její přibližné numerické řešení. Předtím ještě zjednodušeně zavedeme některé pojmy, pomocí nichž se výkony letounů poměřují (jejich určení samozřejmě podléhá jisté standardizaci) • Dolet – vzdálenost, kterou letoun bezpečně uletí • Dostup – maximální výška do které je schopen vystoupat při vodorovném letu • Maximální rychlost – Nejvyšší rychlost, jíž je schopen dosáhnout (většinou se udává v určených letových výškách)

Řešení a simulace letu letadla Vstupní parametry Vzletová hmotnost letadla m0 = 7250 kg, hmotnost letadla při přistání mT = 5250 kg. Nosná plocha S = 23 m2. Tah a spotřeba motoru při jednotlivých režimech: Tah Spotřeba Vypnutý motor 0 kN 0 Μ 100 = 38 kN ∆m100 = −0,4 kg s-1 Plný tah Přídavné spalování

Μ F = 58 kN

∆m F = −3 kg s-1

12

Určení aerodynamických součinitelů Hodnoty bezrozměrných aerodynamických součinitelů naladíme tak, aby určité výkonnostní parametry odpovídaly příslušným reálným hodnotám. Určení aerodynamického součinitele vztlaku: víme, že k „odlepení“ startujícího letounu dojde při rychlosti w = 120 m s-1. V této chvíli tedy bude vztlak větší, než tíhová síla: 1 m0 g < c L ρ 0 S w2 2 odtud plyne pro hodnotu c L 2 m0 g 2 ⋅ 7250⋅ 9,78 cL = =& 0,35 . ; cL = 2 ρ0 S w 1,2255⋅ 23⋅1202 Hodnotu odporového součinitele určíme z maximální rychlosti wMAX = 2125 km h-1, kterou letounu naměřili ve výšce hMAX = 19 km. Nastává silová rovnováha tahu a odporu: 2MF 2 ⋅ 58000 1 2 M F = c D ρ (hMAX ) S wMAX ; c D = =& 0,14 ; cD = 2 2 0,1054⋅ 23⋅ 590,32 ρ (hMAX ) S wMAX

Simulace letu letadla Soustavu rovnic, popisující let letadla budeme řešit přibližnou, numerickou iterativní metodou. (Hodnoty všech pro výpočet potřebných konstant najdeme v předchozím textu. Po celý let ponecháme motor na 100% tahu.). Letoun odstartuje z nulových počátečních podmínek (nultá iterace) s udanou vzletovou hmotností. Budeme iterovat v čase, s konstantním krokem ∆t = 1 s, výpočet zahájíme první iterací. Výpočet bude řízen následující algoritmizací: ALGORITMUS 1: 1. Přechod do další časové hladiny t K = t K −1 + ∆t 2. Hmotnost letounu v K-té iteraci m(t K ) : mK = mK −1 + ∆m100 Tento vztah vyplývá z integrace vztahu (22), resp. úpravě vztahu (23): Pro

∀τ ∈ t K −1 , t K : mK = mK −1 +

tK

∫ ∆m

100

d τ . Provedeme integraci (úbytek hmotnosti je

t K −1

v průběhu časového kroku konstantní: mK = mK −1 + ∆m100 (t K − t K −1 ) = mK −1 + ∆m100 ∆t . Díky volbě časového kroku zjednodušíme tuto rovnici pro výpočet aktuální hmotnosti částečným dosazením: mK = mK −1 + ∆m100 3. Integrací rovnic (26) v časovém kroku zjistíme jednotlivé složky vektoru rychlosti; za předpokladu, že zrychlení je v každém jednotlivém iteračním kroku konstantní: v xK = v xK −1 + a xK −1 a v yK −1 = v yK −1 + a yK −1 Pro

∀τ ∈ t K −1 , t K :

v K = v K −1 +

tK

∫a

K −1

dτ .

Provedeme

integraci

t K −1

v K = v K −1 + a K −1 (t K − t K −1 ) a částečné dosazení: v K = v K −1 + a K −1 4. Integrací rovnic (27) za stejného předpokladu obdržíme vzorce pomocí nichž určíme 1 1 okamžitou polohu letadla: x K = x K −1 + v xK −1 + a xK −1 a hK = hK −1 + v yK −1 + a yK −1 2 2

13

Pro

∀τ ∈ t K −1 , t K :

x K = x K −1 +

tK

∫ (v

K −1

+ a K −1 τ ) d τ .

Provedeme

integraci

t K −1

x K = x K −1 + v K −1 τ + 1 / 2 a K −1 (t K − t K −1 ) a částečné dosazení: x K = x K −1 + v K −1 + 1 / 2 a K −1 2

(

5. Hustota vzduchu v K-té iteraci: ρ K = ρ 0 1 − 0,0965 hK + 0,00325hK *

*2

− 3,7 ⋅10 −5 hK

*3

)

* K −1

nezapomeneme převést dosadit h v jednotkách km! Výpočet realizujeme prostým dosazením známé hodnoty výšky z kroku 4. 1 2 LK = c L ρ K S v xK 6. Aerodynamický odpor a vztlak v K-té iteraci: a 2 1 2 DK = c D ρ K S v xK 2 Výpočet realizujeme prostým dosazením známé hodnoty rychlosti z kroku 3 do Newtonových vzorců. Poznámka.: protože předpokládáme homogenní, rovinné obtékání s relativní rychlostí ve směru osy x, zahrnujeme ve vzorci pouze x – komponentu vektoru rychlosti. Při velkých rychlostech stoupání se tedy dopouštíme chyby!

7. Zrychlení ve směru obou os obdržíme z pohybové rovnice (24): a xK =

M K − DK a mK

− g m K + LK mK 8. Nyní překontrolujeme „kontakt“ letounu se vzletovou drahou: pokud platí, že g mK > LK nahradíme hodnoty y složky rychlosti a polohy nulami: v yK = 0 a hK = 0 . a yK =

9. V případě, že mK ≤ mT motorový let končí, v opačném případě pokračujeme další iterací od bodu 1.. ALGORITMUS 2: 1. Přechod do další časové hladiny t K = t K −1 + ∆t 2. Hmotnost letounu v K-té iteraci m(t K ) : mK = mK −1 + ∆m100 3. Předpokládáme –li, že rychlost stoupání je malá:

(

2

3

)

ρ K = ρ 0 1 − 0,0965 hK −1* + 0,00325hK −1* − 3,7 ⋅10 −5 hK −1* . 1 1 2 2 4. Aerodynamické síly: LK = c L ρ K S v xK −1 a DK = c D ρ K S v xK −1 2 2 M − DK − g mK + LK 5. Zrychlení ve směru obou os z pohybové rovnice: a xK = K , a yK = mK mK 6. Integrací rovnic (26) zjistíme jednotlivé složky vektoru rychlosti; za předpokladu, že zrychlení je v každém jednotlivém iteračním kroku konstantní: v xK = v xK + a xK a v yK = v yK + a yK 1 1 7. x K = x K −1 + v xK + a xK a hK = hK −1 + v yK + a yK 2 2 8. Překontrolujeme „kontakt“ letounu se vzletovou drahou 9. V případě, že mK ≤ mT motorový let končí, v opačném případě pokračujeme další iterací od bodu 1.. Rozdíly mezi hodnotami vypočtené oběma algoritmy jsou minimální.

14

Po provedení numerického výpočtu si necháme vykreslit závislosti rychlosti, výšky a aerodynamických sil na čase: Rychlost 7000 6000 5000 rychlost

4000 vx

3000

vy

2000 1000 0 -1000

0

500

1000

1500

2000

2500

3000

3500

čas

Obr. 11: Závislost jednotlivých komponent rychlosti v x a v y [m s-1] na čase t [s].

výška

Výška 45000 40000 35000 30000 25000 20000 15000 10000 5000 0 -5000 0

500

1000

1500

2000

2500

3000

3500

čas

Obr. 12: Závislost výšky letounu h [m] na čase t [s]. Vztlak a odpor 140000 120000 100000 síla

80000

L

60000

D

40000 20000 0 -20000 0

500

1000

1500

2000

2500

3000

3500

čas

Obr. 13: Závislost aerodynamického vztlaku L a odporu D [N] na čase t [s]. Podíváme –li se na výsledky, zjistíme, že model nefunguje k naší spokojenosti. Ocitli jsme se mimo oblast platnosti některých důležitých předpokladů:

15

•

Velikost rychlosti stoupání je v počáteční fázi letu plně porovnatelná s velikostí dopředné rychlosti letu • Neplatí polynomická závislost hustoty na výšce v celém spektru dosahovaných výšek Navíc dochází k nepředpokládaným oscilacím veličin, které zřejmě způsobuje hmotnost, resp. setrvačnost letounu. Nejvýznamnější chyba modelu ale je nezahrnutí závislosti účinnosti motoru na výšce (všechny veličiny tím rostou nadevšechny meze, odpor je v každém okamžiku menší než tah motoru). Tyto nedostatky se pokusíme odstranit v dalšími – modifikovanými modely.

Modifikovaný model První úprava bude zahrnovat ztrátu účinnosti motoru s rostoucí výškou. Tuto výškovou závislost účinnosti budeme uvažovat nejjednodušší možnou – lineární, procházející body [0m;100% ] a [20000m;50%]. Provedeme numerický výpočet podle stejného schématu, zobrazíme výsledky: Rychlost 300 250

rychlost

200 150

vx

100

vy

50 0 -50

0

500

1000

1500

2000

2500

3000

3500

čas

Obr. 14: Závislost jednotlivých komponent rychlosti v x a v y [m s-1] na čase t [s]. Povšiměme si, že velikost rychlosti stoupání je po vzletu i nadále porovnatelná s dopřednou rychlostí (to nám nevyhovuje). Dobře patrný je okamžik vzletu spojený s nárůstem rychlosti stoupání. A maximální rychlost, která v tomto případě (spíše náhodou, nežli záměrem) dokonce blíží skutečnosti – letoun na plný suchý (100%) tah dosahuje vysokých podzvukových rychlostí.

výška

Výška 16000 14000 12000 10000 8000 6000 4000 2000 0 -2000 0

500

1000

1500

2000

2500

3000

3500

čas

Obr. 15: Závislost výšky letounu h [m] na čase t [s].

16

Vztlak a odpor 140000 120000

síla

100000 80000

L

60000

D

40000 20000 0 0

1000

2000

3000

4000

5000

6000

čas

Obr. 16: Závislost aerodynamického vztlaku L a odporu D [N] na čase t [s]. Ve druhé úpravě se pokusíme reflektovat vliv vysoké rychlosti stoupání. Myšlenka je taková, že rychlost stoupání vyvolává jakousi přídavnou aerodynamickou odporovou sílu T působící proti směru rychlosti stoupání. Rovnici (24a) tedy dostáváme ve tvaru − G + L − T − ma y = 0 Zrychlení ve směru y odtud:

− g m + 1 / 2 c L ρ S v x − 1 / 2 cT ρ S v y 2

ay =

2

; m kde cT je příslušný aerodynamický koeficient. Protože při „kolmém“ pohybu letounu je jeho tvar aerodynamicky nevhodný, odhadneme hodnotu cT = 1. Příslušným způsobem upravíme algoritmus, provedeme numerickou simulaci a zobrazíme výsledky Rychlost 300 250

rychlost

200 150

vx

100

vy

50 0 -50

0

500

1000

1500

2000

2500

3000

3500

čas

Obr. 17: Závislost jednotlivých komponent rychlosti v x a v y [m s-1] na čase t [s].

17

výška

Výška 16000 14000 12000 10000 8000 6000 4000 2000 0 -2000 0

500

1000

1500

2000

2500

3000

3500

čas

Obr. 18: Závislost výšky letounu h [m] na čase t [s]. Vztlak a odpor 140000 120000

síla

100000 80000

L

60000

D

40000 20000 0 0

1000

2000

3000

4000

5000

6000

čas

Obr. 19: Závislost aerodynamického vztlaku L a odporu D [N] na čase t [s]. Jak vidno zavedení přídavného odporu T ovlivňuje pozitivně i fluktuace fyzikálních veličin, přitom ale nemá žádný vliv na dosahované výkony.

Přehled výkonů letadla Výkony jsou stanoveny z výstupních dat numerických simulací druhého, modifikovaného modelu. Poprvé byly data získány při letu na 100% tah, poté se zapnutou forsáží.

plný tah forsáž

doba letu [s] 5 000 670

dolet [km] 1 250 225

dostup [m] 18 300 25 400

maximální rychlost -1 [km h ] 1 095 1 480

18

Závěr Při návrhu matematického modelu letu letadla a jeho následném řešení jsme narazili na jistá úskalí plynoucí buďto z nedostatku dostupných informací, nebo ze složitosti těchto jevů. Tyto kritické situace jsme se snažili vyřešit přiblíženími a zjednodušeními, některé jevy jsme (ne zcela správně) zcela pominuli. Nejvýznamnější vliv na dosahované výsledky má přijetí platnosti Newtonových vzorců (16) pro aerodynamické síly ve všech režimech letu. Tyto vzorce jsou vhodné pro malá Machova čísla. Při rychlostech blížících se rychlosti zvuku a vyšších bychom přesnějších výsledků dosáhli s formulací těchto sil respektující vliv stlačitelnosti a vznik rázových vln. Lepší vyjádření závislosti účinnosti motoru na výšce (tj. závislost reálného tahu motoru na výšce při stálém pracovním režimu motoru) by také zajisté pozitivně ovlivnilo výsledné hodnoty. Stávající výpočet aerodynamických sil (respektující pouze velikost vektoru rychlosti ve směru x) je také problematický, neboť při nenulovém kolmém pohybu již neleží vektor rychlosti nabíhajícího proudu v ose křídla (tětivě profilu). Úhel mezi těmito dvěma entitami nazýváme úhlem náběhu (viz obr.2). V těchto případech se Newtonovy vzorce modifikují tak, že aerodynamické koeficienty jsou závislé na hodnotě úhlu náběhu. Pozn.: V hydromechanice bohužel nelze jednoduše použít princip superpozice a počítat tedy např. s celkovou velikostí vektoru rychlosti v nemodifikovaných Newtonových vztazích s konstantními koeficienty nepovede k lepším výsledkům. Proto je také poměrně diskutabilní zavedení přídavného odporu v „druhé modifikaci“ matematického modelu. Spojitá funkce vázající tahovou sílu se spotřebou paliva v závislosti na režimu chodu motoru by přineslo do modelu mnoho nového. V takto sestaveném modelu by pak bylo možné kupříkladu řídit let letadla plynulým nastavováním tahu motoru tak, aby změny rychlosti stoupání ve stoupavém letu byly co nejmenší. Stanovení aerodynamických součinitelů není jednoznačné. Při výběru jiných výkonových parametrů vedou k různým hodnotám. (Tyto odlišnosti jsou ale provázány s předchozími diskutovanými zjednodušeními). Zajisté zajímavé by bylo získání hodnot koeficientů experimentálně, nebo numerickým výpočtem a coby pevné vstupní hodnoty jejich pomocí určit letounem dosahované výkony. Při řešení složité soustavy rovnic popisující pohyb letadla a síly na něj působící bylo využito primitivní numerické schéma, vnášející do výsledku další chyby. Vhodnější by bylo aplikovat náš model na malý civilní, popř. sportovní letoun létající za bezvětří malými rychlostmi v malých výškách.

19

Použitá literatura Jíra, R., Ticháček, J., Pokorný, V.: Aerodynamika a mechanika letu pro plachtaře. Naše vojsko, SVAZARM, v Praze, 1960. Vykouk, V.: Aerodynamika a mechanika letu pro piloty závěsných kluzáků. ÚV SVAZARMu, v Praze, 1981 Nožička, J., Houbová, R.: Psáno pro L+K: Střípky z aerodynamiky. Časopisy Letectví a kosmonautika č. 7, 8, 9, 10, 11 a 12/2008 Dvořák, R., Kozel, K.: Matematické modelování v aerodynamice. Skriptum ČVUT v Praze, 1996.

20