Teremakusztikai modellezés sztochasztikus módszerrel Bartha Bendegúz VI. vill., Vas Mihály András VI. vill. Tudományos Diákköri dolgozat
Konzulens: Fiala Péter BME Híradástechnikai Tanszék
2008
Tartalomjegyzék 1. Bevezetés
1
2. Sztochasztikus akusztikai térszámítás 2.1. Determinisztikus akusztikai térszámítás . 2.1.1. A hangtér egyenletei . . . . . . . . 2.1.2. Megoldás zárt szobában, módusok 2.1.3. A spektrális végeselem módszer . 2.2. Véletlen rendszerek modellezése . . . . . 2.2.1. A Monte Carlo módszer . . . . . . 2.2.2. A perturbáció módszer . . . . . . . 2.2.3. A módszerek összehasonlítása . . 2.2.4. Egy egyszeru˝ példa . . . . . . . . . 2.2.5. Példa 2 . . . . . . . . . . . . . . . . 2.2.6. Példa 3 . . . . . . . . . . . . . . . . 2.2.7. Konklúzió . . . . . . . . . . . . . . 3. Modellezés 3.1. A spektrális végeselem szoftver . . . . . 3.2. PUS Z TA . . . . . . . . . . . . . . . . . . 3.2.1. Változtatásra váró függvények . 3.2.2. PUS Z TA MCS . . . . . . . . . . . 3.2.3. PUS Z TA perturbáció módszerrel 3.2.4. Módszerek összehasonlítása . . 3.3. Eredmények . . . . . . . . . . . . . . . . 4. A terem mérése 4.1. A mérés el˝okészítése . 4.2. Mérési megfontolások 4.3. A mérés menete . . . . 4.4. A mérés kiértékelése . 4.5. Eredmények . . . . . .
. . . . .
. . . . .
. . . . . i
. . . . .
. . . . .
. . . . .
. . . . .
. . . . .
. . . . .
. . . . .
. . . . . . . . . . . .
. . . . . . . . . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . . . . . . . . .
. . . . . . . . . . . .
3 3 3 5 6 8 9 9 12 12 15 18 19
. . . . . . .
21 21 23 24 24 26 28 29
. . . . .
34 34 35 36 39 40
5. Összefoglalás
42
ii
Ábrák jegyzéke 2.1. Akusztikai térszámítási feladat (forrás: [4]) . . . . . . . . . . 2.2. A (0, 1, 1)-es módus képe . . . . . . . . . . . . . . . . . . . . . 2.3. A véletlen paraméteru˝ rendszerek három lehetséges értelmezése. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.4. A Monte Carlo módszer muködése ˝ . . . . . . . . . . . . . . . 2.5. Tömeg - rugó rendszer: a keresett válasz a tömeg x elmozdulása. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.6. Tömeg - rugó rendszer: véletlen rugómerevség esetén a tömeg kitérésének realizációi (zöld), a perturbáció módszer által szolgáltatott várható érték (kék), a MCS által szolgáltatott várható érték (piros) . . . . . . . . . . . . . . . . . . . . . 2.7. Tömeg - rugó - csillapító rendszer: a keresett válasz a tömeg x elmozdulása. . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.8. Tömeg - rugó - csillapító rendszer: véletlen rugómerevség esetén a tömeg kitérésének realizációi (zöld), a perturbáció módszer által szolgáltatott várható érték (kék), a MCS által szolgáltatott várható érték (piros) . . . . . . . . . . . . . . . . 2.9. Tömeg - rugó - csillapító rendszer: véletlen csillapítás esetén a tömeg kitérésének realizációi (zöld), a perturbáció módszer által szolgáltatott várható érték (kék), a MCS által szolgáltatott várható érték (piros) . . . . . . . . . . . . . . . . . . 3.1. A spektrális végeselem módszer muködését ˝ leíró blokkséma 3.2. Az akusztikai admittancia függése az elnyelési tényez˝ot˝ol . 3.3. A sztochasztikus rendszermodellezés blokksémája Monte Carlo módszer alkalmazása esetén . . . . . . . . . . . . . . . 3.4. A sztochasztikus rendszermodellezés blokksémája perturbáció módszer alkalmazása esetén . . . . . . . . . . . . . . . 3.5. A számítási elrendezés (az IB140 terem felülnézeti képe). Kis pont: vev˝o, nagy pont: vev˝o és gerjesztés . . . . . . . . . 3.6. Az els˝o pár modális koordináta képe. . . . . . . . . . . . . . iii
4 5 8 10 13
15 16
17
20 23 24 25 27 30 31
3.7. A második módushoz tartozó modális koordináta csúcsa. Piros: perturbáció módszerrel a várható érték, kék: Monte Carlo módszerrel a várható érték, zöld: realizációk. . . . . . 3.8. A 21. vev˝opontban a hangnyomáseloszlás. Piros: perturbáció módszerrel a várható érték, kék: Monte Carlo módszerrel a várható érték, zöld: realizációk. . . . . . . . . . . . . . . 3.9. A 21. vev˝opontban három darab módushely. Piros: perturbáció módszerrel a várható érték, kék: Monte Carlo módszerrel a várható érték, zöld: realizációk. . . . . . . . . . . . . 3.10. A 8. vev˝opontban a hangnyomáseloszlás. Piros: perturbáció módszerrel a várható érték, kék: Monte Carlo módszerrel a várható érték, zöld: realizációk. . . . . . . . . . . . . . . . . . 3.11. A 8. vev˝opontban a 37Hz-nél lév˝o módus. Piros: perturbáció módszerrel a várható érték, kék: Monte Carlo módszerrel a várható érték, zöld: realizációk. . . . . . . . . . . . . . . 4.1. 4.2. 4.3. 4.4. 4.5.
Sweepjel id˝otartománybeli képének néhány periódusa . . . Sweepjel spektruma . . . . . . . . . . . . . . . . . . . . . . . . A mérés során használt eszközök . . . . . . . . . . . . . . . . Pillanatkép a mérésb˝ol . . . . . . . . . . . . . . . . . . . . . . A 16-os és 57-es vev˝opontok közti átvitel. Piros: perturbáció módszerrel a várható érték, kék: mérési eredmény. . . . . . . 4.6. A 37-es és 57-es vev˝opontok közti átvitel. Piros: perturbáció módszerrel a várható érték, kék: mérési eredmény. . . . . . . 4.7. A 44-es és 57-es vev˝opontok közti átvitel. Piros: perturbáció módszerrel a várható érték, kék: mérési eredmény. . . . . . .
iv
31 32 32 33 33 36 37 38 38 40 41 41
1. fejezet Bevezetés Az akusztika az ókor óta ismert tudományág, matematikai alapjait, a hanghullámok terjedésének egyenleteit a 19. század óta ismerjük. Ezek analitikus megoldása azonban csak a legegyszerubb, ˝ elméleti jelent˝oségu˝ problémák esetén állítható el˝o, ezért a gyakorlatban csak komoly egyszerusí˝ tések árán használhatóak. Az utóbbi évtizedek számítógépes fejl˝odésének köszönhet˝oen most már lehet˝oségünk van olyan, nagy számításigényu˝ apparátusok használatára, mint a numerikus számítási módszerek. Ezekkel az eszközökkel új távlatok nyithatók, és egy akusztikai problémára pontosabb, a valóságot jobban leíró megoldás adható. A numerikus módszerek a megoldandó egyenletek diszkretizálásával dolgoznak, vagyis az analitikus megoldások véges sok számértékkel jellemzhet˝o közelítéseit, becslését adják. A becslés pontosságának csak a számértékek, az úgynevezett szabadsági fokok száma szab határt, ezért megfelel˝o er˝oforrás birtokában igen pontos számítások végezhet˝oek. Lehet˝ové válik nagy méretu, ˝ bonyolult mechanikai és akusztikai rendszerek modellezés útján való vizsgálata. Problémát okoz azonban, hogy a valóságos problémák numerikus leírásához nagy mennyiségu˝ olyan bemeneti paraméterre (anyagjellemz˝ok, geometria, térbeli pozíciók, stb. . . ) van szükségünk, amit nem feltétlenül ismerünk pontosan. A modellparaméterek bizonytalanságának kétféle oka lehet [2]: 1. A paraméterek determinisztikus jellemzése elvileg lehetséges, de nem ismerjük pontosan az értéküket. Ekkor bizonytalan paraméterekr˝ol beszélünk. Tipikus példa egy tervezés alatt lev˝o autó bels˝o terének akusztikai vizsgálata, ahol az autó méreteit még nem ismerjük pontosan. 2. A paraméterek determinisztikus jellemzése elvileg sem lehetséges. 1
Ezeket változékony vagy véletlen paramétereknek nevezzük. Tipikus példa olyan akusztikai burkolati anyagjellemz˝ok modellezése, melyek akusztikai viselkedése függ a véletlenszeru˝ környezeti behatásoktól (h˝omérséklet, páratartalom), vagy hangterjedési problémák vizsgálata szabadtérben, ahol a hangterjedés paramétereit a véletlenszeru˝ szélmozgás jócskán befolyásolja. Sok esetben a kétféle probléma vegyesen van jelen. A sztochasztikus akusztikai térszámítás feladata az akusztikai hullámok adott geometriai környezetben való terjedésének leírása abban az esetben, mikor a hullámterjedési környezet determinisztikus jellemzése nem lehetséges. A sztochasztikus modellezési eljárások segítségével a bizonytalan bemeneti paraméterek hatása kezelhet˝ové válik, el˝oállíthatók a keresett válasz fontos statisztikai tulajdonságai, és ami talán a legfontosabb, meghatározható a modellezés eredményének bizonytalansága. Jelen dolgozat célja, hogy bemutasson a sztochasztikus akusztikai térszámítási feladat megoldására alkalmas eljárásokat. Minden sztochasztikus módszer megbízható determinisztikus modelleken alapszik. Dolgozatunkban röviden ismertetünk egy determinisztikus akusztikai térszámító módszert és egy azon alapuló szoftvert, majd megmutatjuk, hogy milyen átalakítások segítségével tettük a módszert alkalmassá sztochasztikus paraméterekkel való számításra. Az átalakítások során az eredeti szoftvert „fekete dobozként” kezeltük, ezért muködésé˝ nek részletezése nem tárgya munkánknak. A dolgozat felépítése a továbbiakban a következ˝o: A második fejezet ismerteti a hangtérszámítás feladatát és egy determinisztikus térszámítási módszert, a spektrális végeselem módszert. Ebben a fejezetben foglaljuk össze a sztochasztikus rendszerek modellezésének elméletét, és néhány egyszeru˝ példán keresztül bemutatunk két modellez˝o eljárást, a perturbáció módszert és a Monte Carlo módszert. A harmadik fejezet röviden ismerteti a használt determinisztikus akusztikai térszámító szoftvert, majd bemutatja az abból általunk készített sztochasztikus térszámító programokat. A perturbáció módszeren és a Monte Carlo módszeren alapuló szoftverek részletes ismertetése után tesztek segítségével gyorsasági és pontossági összehasonlítást végzünk. A negyedik fejezet a sztochasztikus modellezés hitelességét igazoló teremakusztikai mérés elméleti és gyakorlati problémáit vizsgálja. Ismertetjük az impulzusválasz mérésének gyakorlati megvalósítását, és egy valós teremben történt mérés eredményeinek segítségével megvizsgáljuk az elméleti modellek és szoftverünk helyességét.
2
2. fejezet Sztochasztikus akusztikai térszámítás 2.1.
Determinisztikus akusztikai térszámítás
Egy forrás által keltett mechanikai rezgést a leveg˝o részecskéi átvéve hanghullám formájában továbbítják. A térnek azt a részét, amelyben a hanghullámok terjednek, hangtérnek nevezzük. A hangtér leírására akusztikai jellemz˝oket használunk. Ezek közül a legfontosabbak a közeg nyomásának pillanatnyi p0 nyugalmi nyomástól vett p(x, t) eltérése, vagyis a hangnyomás és a közeg rezgésben lév˝o részecskéinek v(x, t) sebessége, a részecskesebesség. A hang a szomszédos részecskék közötti energiaátadás útján terjed. A 2.1. ábra alapján látható, hogy a beltéri akusztikai térszámítási feladatok során egy Γa zárt felület által határolt zárt Ωa térfogaton belül, a térfogat egy tetsz˝oleges x pontjában kell meghatároznunk valamely akusztikai jellemz˝o értékét, ill. id˝ofüggését, jelen esetben a lesugárzott p(x) hangnyomást. Ehhez ismernünk kell a peremfeltételeket, vagyis a Γa határoló felület minden pontjában a p(y) hangnyomás vagy a vn (y) normális irányú részecskesebesség értékét (a felület n normálisa a térfogat felé van irányítva).
2.1.1. A hangtér egyenletei Egy forrás által keltett hanghullámok leveg˝oben való terjedését a hangtér egyenletei írják le, megadva kapcsolatot a hangnyomás és a részecskesebesség között a tér egy tetsz˝oleges pontjában. Az egyenletek az alábbi
3
ha vb
ρ0 , c
va n
Ωa Γa
2.1. ábra. Akusztikai térszámítási feladat (forrás: [4]) formában írhatók fel az id˝otartományban [8]: 1 p(x, ˙ t) + ρ0 ∇ · v(x, t) = 0 c2 ρ0 ∇ · v˙ a (x, t) + ∇p(x, t) = 0
(2.1) (2.2)
ahol x a vizsgált tér tetsz˝oleges pontja, t az id˝o, ρ0 a leveg˝o sur ˝ usége, ˝ c pedig a hangsebesség. A változó feletti pont a id˝o szerinti deriválás jele. Ha áttérünk Fourier-transzformációval a frekvenciatartományba, a deriválás iω-val való szorzásra egyszerusödik, ˝ az egyenletek a következ˝o alakra írhatók át: iω ˆ (x, ω) = 0 pˆ(x, ω) + ρ0 ∇ · v c2 ˆ a (x, ω) + ∇ˆ iωρ0 ∇ · v p(x, ω) = 0
(2.3) (2.4)
ˆ (x, ω) pedig a részecsahol pˆ(x, ω) a hangnyomás komplex csúcsértéke, v kesebesség komplex csúcsértéke. A hangtér frekvenciatartománybeli egyenleteib˝ol felírható a Helmholtzegyenlet: ∇2 pˆ(x, ω) + k 2 pˆ(x, ω) = 0 (2.5) ahol k = ω/c az akusztikai hullámszám. Dolgozatban olyan sugárzási problémákkal foglalkozunk, melyekben a rendszer gerjesztése a határoló felület mechanikai rezgése. Amennyiben a határoló felület akusztikai za impedanciája vagy annak reciproka, a ha admittancia ismert, a peremfeltételeket az alábbi módon fogalmazhatjuk meg: pˆ(x, ω) = za (x, ω) [vb (x, ω) − va (x, ω)] n(x) x ∈ Γa
(2.6)
pˆ(x, ω)ha (x, ω) = [vb (x, ω) − va (x, ω)] n(x) x ∈ Γa
(2.7)
4
ˆ b (x, ω) ahol n(x) a térfogat felé mutató normális vektora a Γa felületnek, v jelöli a szerkezeti gerjesztés Fourier-transzformáltját, za (x, ω) és ha (x, ω) pedig a falak akusztikai impedanciájának, illetve admittanciájának frekvencia függ˝o értékei. A (2.6-2.7) egyenletek a szobában kialakuló hangtér és a falak rezgése közti csatolást írják le.
2.1.2. Megoldás zárt szobában, módusok Amennyiben egy zárt szoba válaszát keressük egy forrás által keltett gerjesztésre, a megoldást a hangtér egyenleteinek a Γa zárt felület által határolt, Ωa zárt térfogaton számított megoldása szolgáltatja. Feltételezzük, hogy a szobában létrejöv˝o hangnyomásnak nincs semmilyen hatása az épület rezgéseire. Ez igaz, amennyiben a leveg˝o sur ˝ usége ˝ sokkal kisebb az épületet alkotó anyag sur ˝ uségénél. ˝
2.2. ábra. A (0, 1, 1)-es módus képe Ha a Helmholtz-egyenletet egy Lx × Ly × Lz méretu, ˝ téglatest alakú teremben, merev falat (vn = 0) leíró peremfeltétel mellett oldjuk meg, akkor a terem Ψn (x) módusait kapjuk meg [8], [7]: πlxn x πlyn y πlzn z Ψn (x) = Ψn (x, y, z) = Bn cos cos cos (2.8) Lx Ly Lz Ezek háromdimenziós állóhullámot írnak le, az id˝ofügg˝o eiωt tag nélkül. Ezeket az állóhullámokat hívjuk a terem módusalakjainak, a Bn együttható a módus amplitúdója. A módusokat általában az lxn , lyn , lzn számhármassal szokás megadni, amely a három koordinátatengely irányában kialakuló félhullámok számát jelöli. Egy lehetséges módusalakra mutat példát a 2.2. ábra. A módusok frekvenciája, vagyis a terem sajátfrekvenciái az 5
alábbi módon fejezhet˝ok ki: s 2 2 lny πy 2 c lnx πx lnz πz c flnx lny lnz = + + = 2πkn 2p Lx Ly Lz
(2.9)
A sajátfrekvencia olyan frekvencia, ahol az adott módusalak gerjesztés nélkül kialakulhat és fennmaradhat. Ezeken a frekvenciákon a teremnek, mint rendszernek az átvitelében kiemelések vannak. A terem gerjesztett válaszában az egyes módusalakok a bármilyen frekvencián kialakulhatnak, a módusok szuperpozíciójaként a teremben kialakuló tetsz˝oleges hangnyomáseloszlás felírható. Az efféle számításoknak fels˝o határt szab az úgynevezett fsch Schröder frekvencia [7]: r r c3 T60 fsch = (2.10) 4 ln 10 V Ezen frekvencia felett a módusok száma annyira megszaporodik, hogy nem a terem modális viselkedése dominál, hanem a válasz sokkal inkább közelíthet˝o statisztikusan, nagyfrekvenciás módszerekkel (pl. sugárkövetés).
2.1.3. A spektrális végeselem módszer Akusztikai térszámítási feladatok megoldásának hatékony eszköze a végeselem módszer [11], mely a hangnyomást egyszeru˝ alakfüggvények összegeként írja fel. A spektrális végeselem módszer a módusok összegzésével állítja el˝o a zárt tér egy tetsz˝oleges pontjában a hangnyomás komplex csúcsértékét. pˆ(x, ω) =
∞ X
ˆ n (ω) Ψn (x)Q
(2.11)
n
ˆ n (ω) pedig a módus részeseahol Ψn (x) a zárt tér egy nyomásmódusa, Q dési tényez˝oje (modális koordinátája) az adott frekvencián. Ez a részesedési tényez˝o fejezi ki azt, hogy adott frekvencián az adott módus milyen mértékben vesz részt a hangnyomáskép kialakításában. A szuperpozícióban a módusokat és azok amplitúdóit úgy választjuk meg, hogy azok eleget tegyenek bizonyos ortogonalitási feltételeknek. A módusokra vonatkozó ortogonalitási feltétel: Z Ψn (x)Ψm (x)dΩ = δnm (2.12) Ωa
6
ahol δnm a Kronecker-delta. A módusok térbeli deriváltjaira vonatkozó ortogonalitási feltétel: Z ∇Ψn (x)∇Ψm (x)dΩ = δnm kn km (2.13) Ωa
ahol kn = ωn /c a modális hullámszám. A spektrális végeselem módszer alapegyenlete ezek után a Γa felület által határolt Ωa térfogaton értelmezett Ψn ∇p vektortérb˝ol származtatható, ha alkalmazzuk a Gauss-tételt: Z Z ∇(Ψn (x)∇ˆ p(x, ω))dΩ = − Ψn (x)∇ˆ p(x, ω)n(x)dΓ (2.14) Ωa
Γa
A fenti egyenletek megfelel˝o behelyettesítésekkel és matematikai átalakításokkal (a levezetés megtalálható a [8, 4] forrásokban) az alábbi alakra hozhatók: X ˆ n (ω) + ik ˆ nm (ω) − k 2 Q ˆ n (ω) = kn2 Q Dnm (ω)Q m
ikρ0 c
Z
Ψn (x)ˆ vb (x, ω)n(x)dΓ
(2.15)
Γa
ahol Dnm =
Z
¯ a (x)Ψn (x)Ψm (x)dΓ h
(2.16)
¯ a = ha ρ0 c h
(2.17)
Γa
a csillapításmátrix, és
a normalizált akusztikai admittancia. Ezek után a (2.15) egyenlet tömörebb formája: ˆ =F Λ + ikD − k 2 I Q
(2.18)
ahol Q egy oszlopvektor, amely a modális koordinátákat tartalmazza, Λ a merevségmátrix (diagonális mátrix, f˝oátlójában a modális hullámszámok négyzeteit tartalmazza), D a csillapításmátrix és I az egységmátrix. F a gerjesztés vektor, melynek az n-edik eleme a következ˝oképp írható fel: Z Fn (ω) = ikρ0 c Ψn (x)ˆ vb (x, ω)n(x)dΓ (2.19) Γa
A (2.18) egyenlet segítségével meghatározhatóak a modális koordináták, ha ismertek a módusok és a gerjesztés valamint az impedancia peremfeltételek. A modális koordináták és a módusalakok ismeretében pedig a hangnyomás a tér tetsz˝oleges pontjában számítható. 7
F(ξ(θ))
u(ξ(θ)) K
u(ξ(θ))
F K(ξ(θ))
u(ξ(θ))
F(ξ(θ)) K(ξ(θ))
2.3. ábra. A véletlen paraméteru˝ rendszerek három lehetséges értelmezése.
2.2.
Véletlen rendszerek modellezése
Gyakorlati feladatok megoldása során gyakran el˝oforduló probléma, hogy nincs pontos ismeretünk valamely paraméter értékr˝ol, csak statisztikai tulajdonságai alapján (szórás, várható érték) bizonyos korlátok közé tudjuk szorítani. Ilyenkor a determinisztikus rendszerek modellezésére használt eljárások nem alkalmazhatók, új rendszermodellez˝o módszerekre van szükség: a sztochasztikus vagy véletlen paraméteru˝ rendszerekre. A sztochasztikus rendszer fogalmát több módon is értelmezhetjük, az egyes értelmezési módokat a 2.3. ábra szemlélteti. • Determinisztikus K rendszerátviteli mátrix és véletlen F(ξ(θ)) gerjesztés esetén az u(ξ(θ)) válasz is véletlenszeru˝ lesz. (fels˝o rendszerábra). Itt ξ(θ) egy valószínuségi ˝ változót vagy vektorváltozót jelöl, a véletlenszeruségre ˝ a θ véletlen eseményt˝ol való függés utal. • Determinisztikus gerjesztés és véletlen rendszerátviteli mátrix esetén a válasz szintén véletlenszeru˝ lesz. (középs˝o ábra) • Természetesen véletlen gerjesztés és véletlen rendszerátvitel esetén a válasz ismét csak véletlenszeru˝ lesz. (alsó ábra) Ezek után a sztochasztikus rendszerek modellezésének feladata a következ˝o: Adottak a rendszer K átvitelét és F gerjesztését befolyásoló ξ(θ) véletlen valószínuségi ˝ változók statisztikai tulajdonságai (tipikusan a várható érték és a korrelációmátrix). Ismerjük az átviteli mátrix és a gerjesztés 8
függését a véletlen valószínuségi ˝ változótól. Ezek után feladatunk meghatározni a rendszer F gerjesztésre adott válaszának, az u(ξ(θ)) valószínuségi ˝ vektorváltozónak statisztikai paramétereit.
2.2.1. A Monte Carlo módszer A Monte Carlo módszer [1] olyan, sztochasztikus rendszerek modellezésére használt eljárás, amely a vizsgált rendszer véletlen válaszának mint valószínuségi ˝ változónak statisztikai tulajdonságait (várható értékét, szórását) a „nyers er˝o elvén” határozza meg. Egy adott rendszer esetén el˝oállítja a rendszert és a gerjesztést befolyásoló ξ(θ) véletlen valószínuségi ˝ változó nagy Nreal számú, különböz˝o ξ(θi ) realizációját, majd ezekre különkülön megvizsgálja a rendszer viselkedését, vagyis meghatározza a válasz realizációit. Az eredményekb˝ol, azok statisztikus kiértékelésével meghatározható a válasz valószínuségi ˝ változó várható értéke, szórása, stb. . . A várható érték statisztikai közelítése például: Nreal 1 X ¯ (θ) = E {u(ξ(θ))} ≈ u u(ξ(θi )) Nreal i=1
(2.20)
¯ (θ) átlagérték egy valószínuségi Az u ˝ változó, mely azonos eloszlású valószínuségi ˝ változók átlagaként képz˝odik, így a centrális határeloszlás tétel [9] értelmében √ közel normális eloszlású, szórása pedig az u(ξ(θ)) válasz szórásának Nreal -ed része. Ez az összefüggés fontos útmutató a Monte Carlo módszerben alkalmazandó realizációszám meghatározásában, de a gyakorlati esetekben u(ξ(θ)) szórása nem ismert, így a kell˝o elemszám elérését konvergencviavizsgálattal kell ellen˝orizni [5]. A módszer a realizációk nagy száma miatt rendkívül számításigényes, ezért az eredmények meghatározása gyakran hosszú ideig tart. Egy valószínuségi ˝ modellezéssel megoldandó probléma esetén ezért érdemes megvizsgálni, hogy használható-e másik, kisebb számításigényu˝ módszer a válasz statisztikai tulajdonságainak becslésére.
2.2.2. A perturbáció módszer A perturbáció módszer [5, 10] véletlen rendszerek modellezésére alkalmas módszer, mely abban az esetben alkalmazható el˝onnyel, ha a rendszeregyenlet a következ˝o implicit alakban adott: K(ξ(θ))u(ξ(θ)) = F(ξ(θ)) 9
(2.21)
F(ξ(θ1 ))
F(ξ(θ2 ))
K(ξ(θ1 ))
K(ξ(θ2 ))
u(ξ(θ1 ))
u(ξ(θ2 ))
E {u(ξ(θ))} σ {u(ξ(θ))}
F(ξ(θN ))
K(ξ(θN ))
u(ξ(θN ))
2.4. ábra. A Monte Carlo módszer muködése ˝ ahol K egy négyzetes mátrix. Az egyszeruség ˝ kedvéért feltételeztük, hogy a rendszer átvitele és gerjesztése egyetlen ξ(θ) valószínuségi ˝ változótól függ, melynek várható értéke E {ξ(θ)} = 0. Ezt az általánosság feltételezésének csorbítása nélkül megtehetjük, hiszen a valószínuségi ˝ változót mi választjuk meg. Feltételezzük továbbá, hogy ξ(θ) második momentuma E {ξ 2 (θ)} = 1. Ez a valószínuségi ˝ változó újraskálázásával mindig elérhet˝o. Mivel tudjuk, hogy az F gerjesztés illetve a K átvitel hogyan függ ξ-t˝ol, fel tudjuk írni ezek ξ szerinti 0 körüli Taylor-sorát: 1 ′′ ′ (2.22) F(ξ) = F(0) + Fξ (0)ξ + Fξξ (0)ξ 2 + . . . 2 1 ′′ ′ K(ξ) = K(0) + Kξ (0)ξ + Kξξ (0)ξ 2 + . . . (2.23) 2 Ezekben az egyenletekben ξ hatványainak együtthatói ismertek. Feltételezzük, hogy a válasz is sorba fejthet˝o ξ szerint 0 körül: 1 ′′ ′ (2.24) u(ξ) = u(0) + uξ (0)ξ + uξξ (0)ξ 2 ... 2 Itt ξ együtthatói egyel˝ore ismeretlenek, feladatunk ezek meghatározása. Ha a Taylor-sorokat másodfokú Taylor-polinomokkal közelítjük, akkor a rendszeregyenlet a következ˝oképpen alakul: 1 ′′ 1 ′′ ′ ′ (K(0) + Kξ (0)ξ + Kξξ (0)ξ 2 )(u(0) + uξ (0)ξ + uξξ (0)ξ 2 ) 2 2 1 ′′ ′ = F(0) + Fξ (0)ξ + Fξξ (0)ξ 2 2 10
(2.25)
A megoldást ξ hatványai szerint külön végezzük el. 1. A nulladik hatványhoz tartozó egyenlet: (2.26)
K(0)u(0) = F(0) 2. Az els˝o hatványhoz tartozó egyenlet: ′
′
′
K(0)uξ (0)ξ + Kξ (0)ξu(0) = Fξ (0)ξ
(2.27)
Ezt átrendezve és ξ-vel osztva: ′
′
′
K(0)uξ (0) = Fξ (0) − Kξ (0)u(0)
(2.28)
3. A második hatványhoz tartozó egyenlet: 1 ′′ 1 ′′ 1 ′′ ′ ′ K(0)uξξ (0)ξ 2 ) + Kξ (0)ξuξ (0)ξ + Kξξ (0)ξ 2 u(0) = Fξξ (0)ξ 2 (2.29) 2 2 2 ami átrendezés és 2/ξ 2 -tel való szorzás után: ′′
′′
′
′
′′
K(0)uξξ (0) = Fξξ (0) − 2Kξ (0)uξ (0) − Kξξ (0)u(0)
(2.30)
Az egyenleteket a fenti sorrendben kell megoldani. A megoldás során ′ ′′ el˝oállnak az u(0), uξ (0), uξξ (0) vektorok. Ezek után kiértékelhet˝o az u(ξ(θ)) Taylor-polinomja tetsz˝oleges ξ érték esetén, vagyis ismert a válasz ξ-t˝ol való függésének másodfokú polinomiális közelítése. A válasz várható értéke számítható a nulladfokú és a másodfokú derivált tag és a véletlen változó második momentumából: 1 ′′ ′ 2 E {u(θ)} = E u(0) + uξ (0)ξ + uξξ (0)ξ 2 1 ′′ ′ = u(0) + uξ (0)E {ξ} + uξξ (0)E ξ 2 2 1 ′′ = u(0) + uξξ (0) (2.31) 2 Megfigyelhet˝o, hogy a megoldás átlaga eltér az átlagos paraméteru˝ K(0) rendszer, F(0) átlagos paraméteru˝ gerjesztésre adott u(0) válaszától. Ez az implicit lineáris formában adott rendszerek általános tulajdonsága, annak köszönhet˝o, hogy a mátrixinvertálás nem lineáris muvelet. ˝
11
A válasz korrelációmátrixának közelítése a válasz els˝ofokú Taylor-polinomjának felhasználásával: Ru = E u(θ)uT (θ) n o ′ ′ = E (u(0) + uξ (0)ξ)(u(0) + uξ (0)ξ)T ′
′
= u(0)u(0)T + (u(0)uξ (0)T + uξ (0)u(0)T )E {ξ} ′ ′ +(uξ (0)uξ (0)T )E ξ 2 ′
′
= u(0)u(0)T + uξ (0)uξ (0)T
(2.32)
2.2.3. A módszerek összehasonlítása Látható, hogy a két módszer számításigénye nagymértékben eltér. Míg a Monte Carlo módszer esetén Nreal számú egyenletrendszer megoldása szükséges a válasz várható értékének számításához, a perturbáció módszer ugyanezt a feladatot 3 egyenletrendszer megoldásával elvégzi. Egy adott probléma megoldása során mindenképp ellen˝oriznünk kell, hogy a Monte Carlo szimuláció helyett alkalmazható-e a másik, gyorsabb modellez˝o eljárás. A perturbáció módszer csak olyan esetekben ad helyes eredményt, amikor a paramétereket jellemz˝o ξ valószínuségi ˝ változó és a rendszer válasza között egyszeru˝ összefüggés van, vagyis jó közelítéssel felírható a válasz ξ-t˝ol való függésének alacsony fokszámú polinomiális közelítése. A továbbiakban a módszerek alkalmazását három egyszeru˝ mechanikai példán keresztül illusztráljuk. Az itt kapott eredmények lényegesek lesznek a módszerek akusztikai alkalmazásának szempontjából is.
2.2.4. Egy egyszeru˝ példa Adott a 2.5. ábrán látható egyszeru˝ tömeg-rugó rendszer. A k keménységu˝ rugó egyik vége a merev falhoz van rögzítve, másik fele az m tömeghez kapcsolódik, amelyre egy F determinisztikus er˝o hat. A rugó keménységét a ξ véletlen valószínuségi ˝ változóval jellemezzük, ahol ξ az átlagos rugókeménységt˝ol való véletlenszeru˝ eltérés. k(ξ) = k0 + k1 ξ
(2.33)
Keressük a rendszer válaszát, az m tömeg x elmozdulását. A rendszeregyenlet a következ˝o (az egyszeruség ˝ kedvéért az x válasz ξ-t˝ol való nyilvánvaló függését nem jelöljük): ′′
mx + k(ξ)x = F 12
(2.34)
2.5. ábra. Tömeg - rugó rendszer: a keresett válasz a tömeg x elmozdulása. Ez egy másodrendu, ˝ lineáris, inhomogén, állandó együtthatós differenciálegyenlet, melynek megoldása az id˝otartományban körülményes. Áttérve a komplex frekvenciatartományra, a deriválás s-sel való szorzásra egyszerusödik ˝ és el˝oáll egy sokkal könnyebben kezelhet˝o rendszeregyenlet: mxs2 + k(ξ)x = F (2.35) Kiemelve x-et megkapjuk az implicit rendszeregyenletet: ms2 + k(ξ) x = F ms2 + k0 + k1 ξ x = F
(2.36) (2.37)
ahol ms2 + k0 + k1 ξ = K a perturbáció módszer leírásánál ismertetett átvitel, x a keresett válasz, F pedig a gerjesztés. Ennek megoldása a fent ismertetett lépések alapján: 1. A nulladfokú hatványhoz tartozó egyenlet: x(0) =
F k0 + ms2
(2.38)
mivel ξ = 0. 2. Az els˝ofokú hatványhoz tartozó egyenlet (az átrendezés után): ′
xξ (0) =
−k1 F (k0 + ms2 )2
(2.39)
3. A másodfokú hatványhoz tartozó egyenlet (átrendezés után): ′′
xξξ (0) = 2
13
k12F (k0 + ms2 )3
(2.40)
2.1. táblázat. F m k0 k1 fmax 1 N 0.5 kg 50 N/m 10 N/m 25Hz Ezek után az elmozdulás másodfokú Taylor-polinomos közelítése: ′
′′
x(0) + xξ (0)ξ + xξξ (0)ξ 2 =
−k1 F k12F F + ξ + ξ2 k0 + ms2 (k0 + ms2 )2 (k0 + ms2 )3
A válasz várható értéke: E {x(θ)} = = = =
1 ′′ 2 E x(0) + xξ (0)ξ + xξξ (0)ξ 2 1 ′′ ′ x(0) + xξ (0)E {ξ} + xξξ (0)E ξ 2 2 1 ′′ x(0) + xξξ (0) 2 2 F k12F + E ξ k0 + ms2 (k0 + ms2 )3 ′
(2.41)
Esetünkben ξ legyen -0.5 és 0.5 közötti egyenletes eloszlású véletlen valószínuségi ˝ változó. Az MCS módszer esetén feladatunk el˝oállítani ξ nagyszámú (pl. 20000) realizációját, majd ezekre mind megoldva a rendszeregyenletet, megkapjuk a tömeg kitérésének eloszlását, illetve az eredményeket átlagolva a várható értékét. Perturbáció módszer esetén meg kell határoznunk a gerjesztés, az átviteli mátrix és a válasz deriváltjait((2.38-2.40 egyenletek), majd azok várható értékéb˝ol már számítható a válasz várható értéke a (2.41) képlet szerint. A számításokat a 2.1. táblázatban megadott adatokkal végeztük. A válaszként adódó kitléréseket és a várható érték függvényeket a 2.6. ábra mutatja. Az ábra alapján látható, hogy perturbáció módszer esetén a kitérés spektrumában az r 1 k0 (2.42) f= 2π m rezonanciafrekvencián kiemelés van. Ennek a kiemelésnek a mértéke a különböz˝o ξ értékekt˝ol függ, a rezonanciafrekvencia minden esetben ugyanaz marad. Az MCS módszer esetén a különböz˝o realizációkhoz különböz˝o rezonanciafrekvenciák tartoznak, mivel a rugómerevség szerepel a rendszer sajátfrekvenciájának a képletében. Ezeket átlagolva az eredeti csúcsos rezonancia a várható érték függvényben ellaposodik. Látható, hogy a két 14
kitérés várható értéke[dB]
40 30 20 10 0 −10 −20 −30 1
1.2
1.4 1.6 1.8 frekvencia[Hz]
2
2.2
2.6. ábra. Tömeg - rugó rendszer: véletlen rugómerevség esetén a tömeg kitérésének realizációi (zöld), a perturbáció módszer által szolgáltatott várható érték (kék), a MCS által szolgáltatott várható érték (piros) módszer által meghatározott várható értékek között jelent˝os különbség van. Ez annak köszönhet˝o, hogy a válasz és a valószínuségi ˝ változó közötti összefüggés nem közelíthet˝o alacsony fokszámú polinommal. Ez érthet˝o is, hiszen a rezonanciafrekvencia a rugómerevségt˝ol függ, így az átlagos rendszer rezonanciafrekvenciájának környezetében kiragadott frekvenciaértéken az elmozdulás a rugómerevség nem monoton, gyorsan változó függvénye. Ezért a perturbáció módszer ebben az esetben nem ad helyes eredményt.
2.2.5. Példa 2 Egészítsük ki az el˝oz˝o példát azzal, hogy beiktatunk egy c csillapítót párhuzamosan a rugóval (2.7. ábra). A véletlen valószínuségi ˝ változó most is a k rugókeménység, az els˝o példában ismertetett paraméterekkel. A c csillapítás, az m tömeg és az F er˝o determinisztikus értékek. A rendszert leíró egyenlet a következ˝o: ′′
′
mx + cx + k(ξ)x = F
(2.43)
Ez szintén egy másodrendu, ˝ lineáris, inhomogén, állandó együtthatós differenciálegyenlet, ezért az egyszeruség ˝ kedvéért most is áttérünk a komp15
2.7. ábra. Tömeg - rugó - csillapító rendszer: a keresett válasz a tömeg x elmozdulása. lex frekvenciatartományba. mxs2 + cxs + k(ξ)x = F
(2.44)
kiemelve x-et megkapjuk az implicit alakot: ms2 + cs + k(ξ) x = F
(2.45)
Az MCS módszer esetén feladatunk most el˝oállítani ξ nagyszámú realizációját, majd ezekre mind megoldva a rendszeregyenletet megkapjuk a kitérés eloszlás és a várható értékét. Perturbáció módszer esetén a megoldás során el˝oálló deriváltak: 1. A nulladfokú hatványhoz tartozó egyenlet: F k0 + cs + ms2
x(0) =
(2.46)
mivel ξ = 0. 2. Az els˝ofokú hatványhoz tartozó egyenlet (az átrendezés után): ′
xξ (0) =
−k1 F (k0 + cs + ms2 )2
(2.47)
3. A másodfokú hatványhoz tartozó egyenlet (átrendezés után): ′′
xξξ (0) = 2
k12F (k0 + cs + ms2 )3
(2.48)
A kitérés másodfokú Taylor-polinomos közelítése: ′
′′
x(0) + xξ (0)ξ + xξξ (0)ξ 2 =
F −k1 F + ξ 2 k0 + cs + ms (k0 + cs + ms2 )2 k12F ξ2 (2.49) + (k0 + cs + ms2 )3 16
kitérés várható értéke[dB]
−20 −25 −30 −35 −40 0.5
1
1.5 frekvencia[Hz]
2
2.5
2.8. ábra. Tömeg - rugó - csillapító rendszer: véletlen rugómerevség esetén a tömeg kitérésének realizációi (zöld), a perturbáció módszer által szolgáltatott várható érték (kék), a MCS által szolgáltatott várható érték (piros) F m c k0 k1 fmax 1 N 0.5 kg 1 Ns/m 50 N/m 10 N/m 25Hz 2.2. táblázat. A kitérés várható értéke: 1 ′′ ′ 2 E {x(θ)} = E x(0) + xξ (0)ξ + xξξ (0)ξ 2 1 ′′ ′ = x(0) + xξ (0)E {ξ} + xξξ (0)E ξ 2 2 1 ′′ = x(0) + xξξ (0) 2 2 k12F F + E ξ = k0 + cs + ms2 (k0 + cs + ms2 )3
(2.50)
A számításokat a 2.2 táblázatban látható adatokkal végezve a 2.8. ábrán látható eredményeket kapjuk. A 2.8. ábra alapján látható, hogy a beiktatott csillapítás miatt a rezonanciafrekvencián a várható érték függvények laposabbak lesznek, a kitérések és várható értékük amplitúdója csökken. A két várható érték függvény nem fedi egymást, aminek oka ugyanaz, mint az el˝oz˝o példa esetében: A 17
válasz és a valószínuségi ˝ változó közötti összefüggés a rezonanciafrekvencia környékén most sem közelíthet˝o alacsony fokszámú polinommal.
2.2.6. Példa 3 Legyen most a csillapító c csillapítása a ξ véletlen valószínuségi ˝ változó, amely -0.5 és 0.5 között egyenletes eloszlású, az m tömeg, az F er˝o és a k rugókeménység determinisztikus értékek. Az MCS módszer esetén feladatunk továbbra is el˝oállítani ξ nagyszámú realizációját, majd ezekre mind megoldva a rendszeregyenletet megkapjuk a kitérés eloszlás és a várható értékét. Perturbáció módszer esetén a rendszeregyenlet a következ˝o: ′′ ′ mx + c(ξ)x + kx = F (2.51) ahol (2.52)
c(ξ) = c0 + c1 ξ
Áttérve komplex frekvencia tartományba az implicit rendszeregyenlet a következ˝o: (ms2 + c(ξ)s + k)x = F (2.53) A perturbáció módszer esetén el˝oálló deriváltak: 1. A nulladfokú hatványhoz tartozó egyenlet: F k + c0 s + ms2
x(0) =
(2.54)
mivel ξ = 0. 2. Az els˝ofokú hatványhoz tartozó egyenlet (az átrendezés után): ′
xξ (0) =
−c1 sF (k + c0 s + ms2 )2
(2.55)
3. A másodfokú hatványhoz tartozó egyenlet (átrendezés után): ′′
xξξ (0) = 2
c1 2 s2 F (k + c0 s + ms2 )3
(2.56)
A válasz másodfokú Taylor-polinomos közelítése: ′
′′
x(0) + xξ (0)ξ + xξξ (0)ξ 2 = −c1 sF c1 2 s2 F F + ξ + ξ2 k + c0 s + ms2 (k + c0 s + ms2 )2 (k + c0 s + ms2 )3 18
(2.57)
2.3. táblázat. A 3. példában használt értékek F m c0 c1 k fmax 1 N 0.5 kg 1 Ns/m 0.5 Ns/m 50 N/m 25Hz A kitérés várható értéke: E {x(θ)} = = = =
1 ′′ 2 E x(0) + xξ (0)ξ + xξξ (0)ξ 2 1 ′′ ′ x(0) + xξ (0)E {ξ} + xξξ (0)E ξ 2 2 1 ′′ x(0) + xξξ (0) 2 2 c1 2 s2 F F + E ξ k + c0 s + ms2 (k + c0 s + ms2 )3
′
(2.58)
A számításokat a 2.3 táblázatban található adatokkal végeztük. A választ, vagyis a kitérés értékeket illetve a kitérés várható értékét a 2.9. ábrán láthatjuk. Megfigyelhet˝o hogy a két módszer által szolgáltatott várható érték görbék fedik egymást, a véletlen valószínuségi ˝ változó és a válasz közti összefüggés közelíthet˝o alacsony fokszámú polinommal. Ez érthet˝o, hiszen a csillapítás nagysága az elmozdulásgörbék laposságát vezérli. Az elmozdulás értéke egy adott frekvencián tehát egyszeru, ˝ monoton függvénye a csillapításnak, aminek értelmében a perturbáció módszer jól alkalmazható.
2.2.7. Konklúzió Láttuk, hogy a perturbáció módszer segítségével jól meghatározható véletlen csillapítású rezg˝o rendszerek átvitelének várható értéke. Vegyük észre a 3. példában ismeretett véletlen csillapítású tömeg-rugó-csillapító tendszer és egy véletlen elnyelésu˝ terem közti analógiát. A 2.1.3. fejezetben láttuk, hogy a spektrális végeselem módszer használatakor egy egyszeru˝ geometriájú zárt terembe helyezett forrás hangterének modális koordinátáit a (2.18) szerint számíthatjuk. Ez analóg a (2.53) implicit rendszeregyenlettel, hiszen megtalálható benne a tömegmátrix, a merevségmátrix, a csillapításmátrix, a gerjesztés és a keresett válasz, amely jelen esetben a modális koordináták. Ha feltételezzük tehát, hogy a terem falainak akusztikai impedanciájáról nincs pontos ismeretünk, akkor azt egy véletlen valószínuségi ˝ változóval tudjuk csak leírni. Ennek következtében a (2.16) egyenletben megadott D csillapításmátrix is egy véletlen 19
kitérés várható értéke[dB]
−18 −20 −22 −24 −26 −28 1.2
1.3
1.4
1.5 1.6 frekvencia[Hz]
1.7
1.8
1.9
2.9. ábra. Tömeg - rugó - csillapító rendszer: véletlen csillapítás esetén a tömeg kitérésének realizációi (zöld), a perturbáció módszer által szolgáltatott várható érték (kék), a MCS által szolgáltatott várható érték (piros) valószínuségi ˝ változó lesz. Amint az a 3. példából kiderült, egy ilyen véletlen paraméteru˝ rendszer esetén alkalmazható a perturbáció módszer a válasz, vagyis jelen esetben a modális koordináták meghatározására a lassabb és er˝oforrásigényesebb Monte Carlo módszer helyett.
20
3. fejezet Modellezés A következ˝o fejezet célja, hogy röviden ismertesse a modellezés során átalakított akusztikai térszámító program muködését, ˝ valamint az átalakításának lépéseit, annak érdekében, hogy kezelni tudjon véletlen bemeneti paramétereket. Az átalakítást kétféle módon végeztük el, a már megismert sztochasztikus rendszermodellezési eljárásoknak – Monte Carlo és perturbáció módszer – megfelel˝oen. Az eredmények alapján összehasonlíthatjuk a két módszert, és eldönthetjük, hogy egy adott feladat esetén melyik modell használata el˝onyösebb.
3.1.
A spektrális végeselem szoftver
A P UTA1 névre hallgató, Matlabban írt akusztikai modellez˝o szoftver [3] a spektrális végeselem módszernél leírt számításokat végzi el. A program hatékony kisfrekvenciás teremakusztikai tervezésre ad lehet˝oséget többféle módon is. Ha ismerjük a terem méreteit, a falakat burkoló anyagok admittanciáját és a várható hangforrások elhelyezkedését, a modell segítségével képet kaphatunk a várható hangnyomáseloszlásról, a teremben kialakuló hangképr˝ol. Másik, inverz felhasználási lehet˝oség az, mikor ismeretlen egy terem falainak vagy a falfelületeinek admittanciája, egy akusztikai mérés után azonban rendelkezésünkre áll a terem hangnyomáseloszlása. Ekkor a modellben az ismeretlen admittanciák változtatásával megkaphatunk egy, a mérés eredményéhez nagyon hasonló hangképet, és így becslést adhatunk az ismeretlen admittancia értékekre. Ezek a felhasználási lehet˝oségek feltételezik, hogy a terem bizonyos jellemz˝oir˝ol pontos ismerettel rendelkezünk. Abban az esetben, ha nem 1
Peter’s Universal Toolbox for Acoustics
21
ismerjük pontosan se a hangképet a teremben, se a burkolatok anyagjellemz˝oit, a numerikus teremakusztikai modellez˝o eljárások hatékony segítséget nyújthatnak. A szoftvert kiegészítve sztochasztikus rendszermodellezési funkciókkal, egy olyan eszköz kerül a kezünkbe, amellyel kisfrekvencián, az egyszeru˝ geometriájú termekben megbecsülhetjük a hangnyomás várható értékét tetsz˝oleges pontban, ha a terem falainak elnyelése nem ismert pontosan, csak statisztikus jellemz˝oi által (szórás, várható érték) bizonyos korlátok közé szorítható. Ez által egy, a gyakorlati tervezés során is jól alkalmazható programot kapunk. A következ˝o fejezetben részletezett átalakítások megértéséhez feltétlenül szükséges egy rövid ismertet˝o a program által használt függvények muködésér˝ ˝ ol. A program blokksémája a 3.1. ábrán látható. 1. Módusok meghatározása (modes függvény): A terem Li méretei, a hangsebesség és a fels˝o határfrekvencia alapján meghatározza a terem módusalakjait és azok sajátfrekvenciáit. 2. Csillapításmátrix meghatározása (damping függvény): A fal ha admittanciáinak megadása után meghatározza falfelületek za akusztikai impedanciáját, majd a D csillapítási mátrixot. Az elnyelés lehet térben homogén, négyszögletes területeken belül homogén illetve tetsz˝oleges eloszlású. 3. A gerjesztés meghatározása (excitation függvény): A gerjesztés meghatározása során meg kell adnunk annak pozícióját, valamint hogy milyen jellegu˝ a forrásunk. Lehet felület-, illetve pontforrás, amelyet elhelyezhetünk a falakon vagy a bels˝o térfogatban valahol. A program képes kezelni egyszerre több különböz˝o forrást is. 4. A modális koordináták meghatározása (solve függvény): Meghaˆ n modális koordinátákat a leveg˝o hullámimpedanciája, a tározza Q csillapítási mátrix, a frekvenciavektor, és az excitation függvény eredményeként el˝oálló gerjesztési vektor ismeretében. Tehát a solve függvény oldja meg a (2.18) egyenletben bevezetett ˆ =F Λ + ikD − k 2 I Q egyenletrendszert.
5. A hangnyomás meghatározása (response függvény): Az xr vev˝opontok megadása után a modális koordináták és a módusalakok alapján kiszámolja a hangnyomás komplex csúcsértékét a vev˝opontokban a kívánt frekvenciákon. 22
ˆ (x, ω) v Lx , Ly , Lz c modes
fmax
Ψn (x), kn
excitation
F(ω)
D
damping
solve
ha ˆ n (ω, θi ) Q
response xr p ˆ (xr , ω)
3.1. ábra. A spektrális végeselem módszer muködését ˝ leíró blokkséma
3.2.
PUS Z TA
Ebben a fejezetben bemutatjuk, hogy a már megismert spektrális végeselem modellez˝o program, mely részeit és hogyan kell átalakítani ahhoz, hogy a támasztott követelményeknek megfeleljen, vagyis az admittanciával mint véletlen valószínuségi ˝ változóval képes legyen dolgozni. Az így implementált új programot PUS Z TA névre kereszteltük, melyben az „Sz” a sztochasztikus üzemmódra utal. A megvalósítás során a már megismert két módszerrel, a Monte Carlo és a perturbáció módszerrel fogunk dolgozni. A valóságot jól leíró módon, valószínuségi ˝ változónak az α elnyelési tényez˝ot, illetve az abból a 3.1 képlettel számolható akusztikai admittanciát választjuk. Az erlnyelési tényez˝o gyakran használt anyagjellemz˝o menyiség a teremakusztikában. Értéke a falra mer˝olegesen bees˝o síkhullám energiájának és a fal által elnyelt hangenergiának arányát adja meg. A 3.2 23
ábrán láthatjuk, hogy az elnyelési tényez˝o és az akusztzikai admittancia összefüggése közel lineáris. Az egyszerusítés ˝ érdekében az elnyelési tényez˝o térben homogén a felületen. 0.03
Relatív admittancia[−]
0.025
0.02
0.015
0.01
0.005
0 0
0.02
0.04 0.06 Elnyelési tényezõ[−]
0.08
0.1
3.2. ábra. Az akusztikai admittancia függése az elnyelési tényez˝ot˝ol √ 1 1− 1−α √ ha = ρ0 c 1 + 1 − α
(3.1)
3.2.1. Változtatásra váró függvények A kétféle módszer különböz˝o megvalósítást kíván, de mindkett˝oben közös, hogy a f˝oprogram függvényhívásai, illetve paramétermeghatározሠn (ω) modális koorsai, a D csillapítás mátrixot meghatározó damping, a Q dinátákat meghatározó solve, valamint a választ, vagyis a hangnyomás pˆ(x, ω) komplex csúcsértékét számító response függvények a változásokban érintett részek.
3.2.2.
PUS Z TA MCS
A Monte Carlo módszert használva, a számításaink pontosságának csupán a realizációk száma szab határt. Így az MCS módszerrel dolgozó programot nyugodtan tekinthetjük a referenciánknak, mellyel megállapíthatjuk a perturbáció módszer használatának jogosságát. 24
ˆ (ω) v Lx , Ly , Lz c fmax
modes
Ψn (x), kn
excitation
F(ω)
damping
D
d
ha
ˆ n (ω, θi ) Q
response xr pˆ(xr , ω, θi )
3.3. ábra. A sztochasztikus rendszermodellezés blokksémája Monte Carlo módszer alkalmazása esetén A referencia szerepért azonban áldozatokat kell hoznunk. Mert igaz bár, hogy a pontosságnak elméleti korlátja nincs, a gyakorlatban azonban komoly problémát okoz a futási id˝o, melyet a realizációk nagy száma okoz. A realizációk számának megválasztása nem egyszeru˝ feladat, hiszen a túl kevés realizáció hamis eredményhez vezet, míg a kelleténél több realizáció feleslegesen nagy futási id˝ovel jár. Annak a megállapítására, hogy mekkora a szükséges realizációszám, jelen esetben egy igen egyszeru˝ módszert alkalmaztunk. A realizációk számának növelésével a vizsgált értékek, jelen esetben a várható érték és a szórás, egy értékhez, az elméleti eredményhez konvergál. Iterációval megkerestük azt a legkisebb értéket, ahonnan a realizációk számának növelése már nem okozott változást a vizsgált értékekben. A kapott eredményeket a továbbiakban helyesnek fogadjuk el, és a kés˝obbiekben ezekkel hasonlítjuk össze a perturbáció módszerrel kapott eredményeinket. A konkrét programozási feladat – a kód megfelel˝o megváltoztatása 25
– ebben az esetben nagyon egyszeru, ˝ a függvényeket nem kell módosítani. A változtatás tulajdonképpen abból áll, hogy a realizációk számának megfelel˝oen, a függvények meghívását ciklusba kell szervezni, és minden egyes realizációra lefuttatni. A változások az elnyelési tényez˝o meghatározásánál kezd˝odnek, ahol az egyszeri számbevitel helyett nagyszámú realizáció vektroba rendezése történik. A következ˝o változás a D csillapításmátrix meghatározását végz˝o damping függvény meghívásánál történik. A csillapításmátrixot meghatározó egyenletben, ha a za akusztikai impedancia, illetve annak reciproka, a h admittancia a felületen homogén, vagyis nem függ x-t˝ol, akkor kiemelhet˝o az integrál elé. Z ¯ Dnm = ha (x) Ψn (x)Ψm (x)dΓ (3.2) Γa
Így nagyban leegyszerusíthetjük ˝ a dolgunkat, hiszen a nagyszámú realizáció mellé nem szükséges nagyszámú csillapításmátrix legenerálása és tárolása, hiszen az admittanciával való szorzást a kés˝obbiekben is elvégezhetjük, ezt kihasználva a damping függvény meghívásakor az admittanˆ n modális kordináták számícia értékét egyel˝ore 1-re állítjuk. Innent˝ol a Q tásáig nem történik változás. A modális koordináták számításánál egy for ciklust hívunk meg, ahol minden egyes realizációra meghatározzuk a modális koordinátákat a következ˝o módon: el˝oször legeneráljuk az éppen aktuális csillapításmátrixot, majd ezt helyettesítjük be a változatlanul hagyott solve függvénybe. E módon megkaptuk a modális koordináták sokaságát minden realizációra, ezt egy háromdimenziós mátrixban tároljuk el. Ennek a háromdimenziós mátrixnak egy-egy síkja jelent egy realizációt. A megoldás meghatározásánál a response függvényt szintén egy for ciklusba foglaljuk, és minden egyes síkra elvégezzük, így megkaptuk a válaszok sokaságát. Látható, hogy a Monte Carlo módszer sok számítással jár, miközben nagy adathalmazokat kell eltárolnunk. Ennek megfelel˝oen futási ideje rendkívül hosszú is lehet. El˝onye viszont, hogy egyszeruen ˝ lekódolható.
3.2.3.
PUS Z TA perturbáció módszerrel
Mint azt a 2.2. fejezetben említettük, a perturbáció módszer alkalmazhatóságának határt szab, hogy milyen a válasz függése a véletlen valószínu˝ ségi változótól. A szimuláció során a modális koordinátákat fogjuk a Taylor-polinomjukkal közelíteni. A kés˝obbiekben meg fogjuk mutatni, hogy ebben az esetben a módszer alkalmazása lehetséges.
26
vˆ(x, ω) Lx , Ly , Lz Ψn (x), kn
c modes
fmax
excitation
F(ω)
damping
E[q]
D
E[p] solve
pˆ(xr , ω, θi )
2
σ [ˆ q (ω)] 2 σ [p]
ha
response xr
3.4. ábra. A sztochasztikus rendszermodellezés blokksémája perturbáció módszer alkalmazása esetén A megvalósítást itt is az elnyelési tényez˝ok sokaságának el˝oállításával kezdjük, ennek azonban itt az a célja, hogy numerikusan meghatározzuk a ξ valószínuségi ˝ változó E {ξ} várható értékét és E {ξ 2 } második momentumát. Mivel ez egy egyszeru˝ muvelet, ˝ így a kelleténél jóval nagyobb realizációszámnál is rendkívül gyorsan kiszámolható. Ezt követ˝oen a h(ξ) akusztikai admittancia, és annak ξ szerinti els˝o és másodfokú deriváltjainak meghatározása következik a ξ = 0-helyen. A damping függvényt az MCS módszernél megismert okokból az ott ismertetett módon hívjuk meg, és ahhoz hasonlóan a modális koordináták meghatározását végz˝o részig innent˝ol fogva nincs változás. A perturbáció módszernél azonban nem alkalmazunk ciklusokat, hanem két függvényt is módosítunk a számítás elvégzésének érdekében. El˝oször is megfordítjuk a sorrendet, és el˝oször a response függvényt hívjuk meg, ennek azonban csökkentjük a funkcióit, mert nem a nyomás komplex csúcsértékét határozza meg, hanem csupán a Ψn (x) nyomásmódusokat a kiválasztott x pozíciókban. Ezeket bemenetként meg kell ad27
nunk a solve függvénynek, mely teljes mértékben f˝oszerepl˝ové lép el˝o, ˆ n modális koordináták, és a pˆ nyomás komplex hiszen meghatározza a Q csúcsértékének várható értékét és szórását. A függvény a frekvenciafelbontásnak megfelel˝oen egy for cikluson belül polinomiális közelítéssel számol várható értéket és korrelációmátrixot, melynek segítségével a szórás számolható. A válaszszámítás hozzáadása ehhez a függvényhez azért volt indokolt, mivel a korrelációmátrixok tárolása feleslegesen megbonyolította volna a feladatot. A (3.3) egyenletben látható, hogy a korrelációmátrix diagonálisában a valószínuségi ˝ változó második momentuma található [9], ami behelyettesítve a (3.4) egyenletbe, számolható a szórásnégyzet. A perturbáció módszer tárgyalásánál láthattuk, hogy a (2.32) segítségével számolhatjuk a modális koordináták korrelációmátrixát. A nyomásválasz korrelációmátrixát a a (3.3) egyenlet alapján a (3.5)-ben látható módon száP ˆ míthatjuk. Mivel pˆ(x, ω) = ∞ n Ψn (x)Qn (ω), és Ψ konstans, ezért a várható érték képzésb˝ol kiemelhet˝o, így a (3.6)-ben látható módon megkapjuk a válasz korrelációmátrixát, és abból már számítható a válasz szórása. Rξ = E ξ(θ)ξ T (θ) E {ξ12 (θ)} E {ξ1 (θ)ξ2 (θ)} E {ξ2 (θ)ξ1 (θ)} E {ξ22 (θ)} = .. .. . . E {ξn (θ)ξ1 (θ)} E {ξn (θ)ξ2 (θ)}
. . . E {ξ1 (θ)ξn (θ)} . . . E {ξ2 (θ)ξn (θ)} (3.3) .. ... . ...
E {ξn2 (θ)}
D {ξ(θ)}2 = E ξ 2 (θ) − m2ξ
(3.4)
Rp = ΨE Q · QT ΨT = ΨRQ ΨT
(3.6)
Rp = E p · p T
(3.5)
3.2.4. Módszerek összehasonlítása
Megvizsgálva a módszerek adta lehet˝oségeket, megállapíthatjuk, hogy a perturbáció módszer drasztikusan kisebb futási id˝oket eredményez, mint az MCS módszer. Azonban ahhoz, hogy ezt az el˝onyt kiélvezhessük, bonyolultabb programkóddal és a megvalósítás el˝ott komoly megfontolásokkal kell számolnunk. Azonban meg fogjuk mutatni, hogy ebben a konkrét esetben a perturbáció módszer használata lehetséges, s˝ot ajánlott, hiszen vele egy sokkal rugalmasabban, és könnyebben használható eszközt kapunk. 28
3.1. táblázat. Az I.B.140-es terem méretei Hosszabbik oldal Rövidebbik oldal Magasság 7.66m 5.84m 2.93m ¯ a véletlen valószínuségi 3.2. táblázat. A h ˝ változó paraméterei ¯a ¯a E h σ h 0.0077 0.0045
3.3.
Eredmények
Vizsgáljuk most meg a a modellezési eljárásokkal kapott eredményeinket. A modell egy valós tanteremr˝ol, a BME I épületének I.B.140-es termér˝ol készült, így mérésekkel is igazolhatjuk a modell helyességét (lásd a 4. fejezetet). A terem paraméterei közül a mérete, és felületein az elnyelési tényez˝onek az ismerete szükséges. Bár a teremben többféle felület található, az egyszerubb ˝ számítás érdekében homogén elnyelésu˝ felületetként kezeljük a falakat. Véletlen változónak a ha akusztikai admittanciát választjuk, ami könnyedén megfeleltethet˝o az elnyelési tényez˝onek, ahogy azt a 3.2. fejezetben láttuk. ¯ a = ρ0 cha relatív admittanA modellezés során ha -t – pontosabban a h ciát egyenletes eloszlású véletlen változónak definiáltuk, a 3.2 táblázatban található paraméterek szerint. A számítási elrendezés a 3.5. ábrán látható. Gerjesztésnek egy pontforrást alkalmaztunk, mely a terem egyik sarka közelében helyezkedik el. A nyomásválaszt a terem 9 × 7-es hálón felvett 63 különböz˝o pontjában számítottuk. Ezek után vizsgáljuk meg a modellezési eredményeinket! A 3.6. ábrán láthatjuk a modális koordinátákat, melyek szemléletesen mutatják, hogy az egyes módusok adott frekvencián mekkora részesedést vállalnak a válasz kialakításában. Látszik, hogy minden egyes modális koordináta felfogható egy egyszeru˝ mechanikai rezg˝orendszer válaszaként. A 3.7 ábrán a második modális koordináta csúcsát látjuk kinagyítva. Ábrázoltuk a Monte Carlo módszerrel kapott realizációkat, azok várható értékét és a perturbáci módszer által adott várható érték görbét. Id˝ohiány miatt a Monte Carlo módszert nem tudtuk kell˝o finomságu frekvenciaskálával futtatni. Megfigyelhetjük azonban, hogy a két várhatóérték görbe együtt fut, tehát a perturbáció módszer alkalmas ebben az esetben a problémamegoldásra. A 3.8 és a 3.9 ábrákon a 21-es számú, vagyis az x = 1.9 m, y = 5.74 m 29
6 5
y [m]
4 3 2 1 0 0
2
4 x [m]
6
8
3.5. ábra. A számítási elrendezés (az IB140 terem felülnézeti képe). Kis pont: vev˝o, nagy pont: vev˝o és gerjesztés vev˝opontban hogyan alakul a hangnyomáskép, illetve a kritikus pontokon, a módushelyeken hogyan viselkednek a modellez˝o eljárások. Itt is láthatjuk, hogy a két modellez˝o eljárás várható értéke együtt fut, mint ahogyan ugyan ezt figyelhetjük meg a 8-as mér˝oponthoz (x = 1 m, y = 0.2 m) tartozó azonos jellegu˝ a 3.10 és a 3.11 ábrákon is. Ezek után bátran kijelenthetjük, hogy adott modellezési körülmények között a perturbáció módszer maradéktalanul alkalmazható. A gyakorlati tapasztalatok azt mutatják, hogy azonos feltételek mellett a Monte Carlo szimulációnak több órányi a futásideje, míg a perturbáció módszerét alkalmazó eljárás egy-két perc alatt a Monte-Carlo szimulációnál nagyobb frekvenciafelbontásban képes ugyanazt az eredményt szolgáltatni. A következ˝okben azt vizsgáljuk, hogy ez a modellez˝o eljárás mennyire fedi a valóságot, vagyis mennyire alkalmazható a gyakorlatban.
30
80 70 60
Amplitúdó[dB]
50 40 30 20 10 0 −10 −20 0
20
40
60 80 Frekvencia[Hz]
100
120
140
3.6. ábra. Az els˝o pár modális koordináta képe.
76 75 74
Amplitúdó[Hz]
73 72 71 70 69 68 67 22
22.2
22.4 22.6 Frekvencia[Hz]
22.8
23
3.7. ábra. A második módushoz tartozó modális koordináta csúcsa. Piros: perturbáció módszerrel a várható érték, kék: Monte Carlo módszerrel a várható érték, zöld: realizációk. 31
70
60
Amplitúdó[dB]
50
40
30
20
10
0 0
20
40
60 80 Frekvencia[Hz]
100
120
3.8. ábra. A 21. vev˝opontban a hangnyomáseloszlás. Piros: perturbáció módszerrel a várható érték, kék: Monte Carlo módszerrel a várható érték, zöld: realizációk.
58
56
Amplitúdó[dB]
54
52
50
48
46
56
58
60
62 64 Frekvencia[Hz]
66
68
70
3.9. ábra. A 21. vev˝opontban három darab módushely. Piros: perturbáció módszerrel a várható érték, kék: Monte Carlo módszerrel a várható érték, zöld: realizációk. 32
70
60
Amplitdó[dB]
50
40
30
20
10
0 0
20
40
60 80 Frekvencia[Hz]
100
120
3.10. ábra. A 8. vev˝opontban a hangnyomáseloszlás. Piros: perturbáció módszerrel a várható érték, kék: Monte Carlo módszerrel a várható érték, zöld: realizációk.
62
Amplitdó[dB]
60
58
56
54
52
50
36.2
36.4
36.6
36.8 37 37.2 Frekvencia[Hz]
37.4
37.6
37.8
3.11. ábra. A 8. vev˝opontban a 37Hz-nél lév˝o módus. Piros: perturbáció módszerrel a várható érték, kék: Monte Carlo módszerrel a várható érték, zöld: realizációk. 33
4. fejezet A terem mérése Annak érdekében, hogy a modellezési eljárásunkat ellen˝orizhessük, a modellezett teremben méréseket végeztünk. A mérés során ugyanazt az elrendezést használtuk, mint amit az el˝oz˝o fejezetben leírt számítási példában: 63 vev˝opontot és egy gerjesztést helyeztünk el a tanteremben.
4.1.
A mérés elokészítése ˝
Mivel a modellt a lehet˝o legnagyobb mértékben szeretnénk leegyszerusí˝ teni, ezért a mérést is üres teremben kellett végrehajtani. A vev˝opontok nagy száma miatt több mikrofon párhuzamos használata indokolt, így a konkrét esetben hét vev˝opontot határoztunk meg a terem rövidebbik oldala mentén, ez megfelel egy mikrofonsornak. A mikrofonsor jeleit a terem hosszabbik oldala mentén kilenc helyen rögzítettük, így a teremben egyenletesen elosztottuk a hatvanhárom vev˝opontot. A mikrofonok a magassági tengelyb˝ol egy síkot metszettek ki, mely párhuzamosan helyezkedett el a padlóval. A gerjesztést egy hangsugárzó szolgáltatta, ezt az egyik sarokban, egy vev˝opont alatt helyeztük el. A mérési elrendezést az el˝oz˝o fejezetben ismertetett a 3.5. ábra mutatja. A mérés el˝okészítésekor fontos szempont az id˝opont megválasztása. Szabad téri mérésnél figyelembe kell venni a légköri változások, és ezáltal az akusztikai paraméterek változását. Jelen esetben beltéri mérésr˝ol van szó, ezért ett˝ol a hatástól eltekinthetünk, azonban egy másik fontos jelenségr˝ol nem. A termünk egy nagy, sok ember által használt, sur ˝ u˝ forgalmú útak által határolt épületben van, ezért a háttérzajokról nem szabad megfeledkeznünk. Ezek a háttérzajok nagyban ronthatják a mérésünk eredményeit, mivel a megemelkedett alapzaj, vagy másképpen kifejezve a folyamatosan jelenlév˝o haszontalan háttérenergia rontja a mérés jel-zaj vi34
szonyát. Ezek a háttérzajok származhatnak emberek által keltett zajokból – ilyen például a beszéd, a mozgás –, az épület gépészeti zajaiból – ez lehet például klímaberendezés, lift –, valamint az épületbe beszur˝ ˝ od˝o küls˝o zajokból – ez jellemz˝oen forgalmi zajt jelent. Láthatjuk, hogy ezen zajok többsége a nappali általános munkavégzés eredménye, ezért célszeru˝ a mérésnek éjszakai id˝opontot választani. Ezzel megkíméljük magunkat a háttérzajok kellemetlen hatásától, hiszen az épület jellemz˝oen üres, így akár a gépészet is kikapcsolható, a forgalom pedig jóval a nappali szint alatt van. Ezen kívül, bár ez nem közvetlenül kapcsolódik a méréshez, nem zavarjuk az épületben dolgozókat, hiszen nem feledkezhetünk meg arról, hogy a mérés során nagy hangnyomásszintu˝ gerjeszt˝ojelekkel dolgozunk. Így ebben a konkrét esetben is az éjszakai mérés mellett döntöttünk.
4.2.
Mérési megfontolások
A mérés során impulzusválaszokat szeretnénk mérni, melyek segítségével átviteli függvényeket határozhatunk meg. Ezeket az átviteli függvényeket hasonlíthatjuk össze a modell segítségével számolt átviteli függvényekkel. Az impulzusválasz meghatározható közvetlen és közvetett módszerekkel [6] , [7]. A közvetlen módszer azt jelenti, hogy megpróbálunk impulzusgerjesztést el˝oállítani, és erre mérjük a választ, így a mért jel az impulzusválasz lesz. Ezeknek a módszereknek két kritikus pontja van, az egyik maga a jel el˝oállítása, a másik pedig a mérés jel-zaj viszonya. A jel el˝oállítása többféle módon történhet, például pisztolylövés, szikrakisülés, hangsugárzóval történ˝o impulzuskijátszás, de valójában egyik módszerrel sem érhetünk el megfelel˝o spektrumú, tehát a kell˝o frekvenciasávban egyenletes átvitelu˝ gerjeszt˝o jelet. A megfelel˝o jel-zaj viszony elérése is problémás, hiszen ezekkel a gerjeszt˝o jelekkel nem tudunk megfelel˝o mennyiségu˝ energiát bevinni a rendszerbe. A közvetett módszerek lényege abban rejlik, hogy nem impulzusgerjesztéssel dolgozunk, hanem választunk egy, az igényeinknek megfelel˝o jelet, majd az arra kapott választ különböz˝o matematikai apparátusok segítségével impulzusválasszá alakítjuk. A közvetett módszereket két csoportba sorolhatjuk, a szélessávú, vagy keskenysávú módszerekhez. A szélessávú módszer lényege abban rejlik, hogy egy olyan jelet választunk gerjeszt˝o jelként, mely a számunkra érdekes sávban minden frekvenciakomponenst azonos energiával tartalmaz. Jellemz˝oen ez a gerjeszt˝o jel a fehérzaj. A szélessávú módszer el˝onye, hogy nagyon gyors mérést tesz lehet˝ové, hátránya, hogy a jel-zaj viszonya alacsony, hiszen az egyes frekvenicakomponensek egyenként kevés energiát közvetítenek. 35
A keskenysávú módszereknél egy vagy kevés frekvenciakomponenst tartalmazó jelet használunk, melynek a frekvenciáját változtatva végigmegyünk a kívánt frekvenciatartományon. A módszer el˝onye, hogy tulajdonképpen bármekkora energiát közvetíthetünk az egyes frekvenciákon, így a jel-zaj viszonyunk jó lesz, hátránya azonban, hogy a mérés hosszú ideig tart. Az IB140-es terem mérésénél, a minél pontosabb eredmény elérésének érdekében közvetett keskenysávú mérést végeztünk, a gerjeszt˝o jelünk egy exponenciális sweepjel volt, annak érdekében, hogy a kritikus alacsony frekvenciákon a lehet˝o legtöbb energiát vigyük a rendszerbe. A jel id˝ofüggvényét [7] a (4.1) képletben láthatjuk: ! i ωalso · T h T1 log( ωωfelso ) also sweep = A · sin −1 (4.1) · e felso log( ωωalso ) ahol A az amplitúdó, ωalso az alsó határfrekvencia, ωalso a fels˝o határfrekvencia, T pedig a sweep jel hossza. A mérés kiértékelésénél – lásd a 4.4. fejezetet – át kell térnünk a frekvenciatartományba, ehhez pedig a jelünket ablakoznunk kell, ennek hatását láthatjuk a 4.1 ábrán. 1
Relatív amplitúdó
0.5
0
−0.5
−1 0
0.02
0.04
0.06
0.08
0.1 Idô[s]
0.12
0.14
0.16
0.18
0.2
4.1. ábra. Sweepjel id˝otartománybeli képének néhány periódusa
4.3.
A mérés menete
A mérés során a gerjeszt˝o jel hossza a bevitt energiával, így a jel-zaj viszonnyal arányos. Annak érdekében, hogy a váratlanul fellép˝o véletlen zajokat is kiszurjük ˝ a mérésb˝ol, átlagolást alkalmazhatunk oly módon, hogy 36
0
Amplitúdó[dB]
−10
−20
−30
−40
−50 0
100
200
300 Frekvencia[Hz]
400
500
600
4.2. ábra. Sweepjel spektruma a gerjeszt˝o jelet többször megismételjük, így a rá adott választ többször rögzíthetjük, melyeket átlagolva, a változó véletlen zajok hatása nem jelenik meg az eredményünkben. Két gerjeszt˝o jel között eltelt id˝o [6] határozza meg a mért impulzusválasz hasznos hosszát, vagyis a mért impulzusválasz olyan hoszzú lesz, amekorra két gerjeszt˝o jel megszólalása között eltelt id˝o. A fentiek ismeretében a mérés során négyszeres ismtélésu, ˝ hatvan másodperc hosszú, az ismétlések között három másodperc szünetet tartalmazó sweep jelsorozattal dolgoztunk. Mivel a mért terem egy átlagos tanterem, a három másodperces impulzusválasz hossz nagy biztonsággal elegend˝o. Miután a mérés el˝okészítése megtörtént, és felállításra került az els˝o mikrofonsor, megkezd˝othetett a konkrét mérés. Miután hosszú sweepjelekkel dolgoztunk, négyszeres ismétl˝odéssel, és a mikrofonsorokat kilenc helyre pozícionáltuk, ezért az effektív mérés id˝oben szintén hosszú. Azonban, mivel bels˝o térben dolgoztunk, azonos körülményeket feltételezhetünk a mérés teljes hossza alatt, a küls˝o paraméterek nem változtak olyan mértékben, hogy az az eredményeinket meghamisítaná. Mér˝orendszerünk kondenzátormikrofonokból, az azokhoz alkalmazandó töltéser˝osít˝ot magában foglaló mér˝oer˝osít˝ob˝ol és egy hangkártyából állt. A mér˝oeszközök pontosabb jellemzése a 4.1. táblázatban található meg. A mér˝oeszközöket és a mérés során készült pillanatfelvételt a 4.3-4.4. ábrákon mutatjuk meg.
37
4.3. ábra. A mérés során használt eszközök
4.4. ábra. Pillanatkép a mérésb˝ol
38
4.1. táblázat. A mérés során használt eszközök Gyártó Motu
Típus 828mkII
PCB
482A20
PCB
130P10
Genelec Genelec Dell
7070A
Eszköz 8 csatornás küls˝o hangkártya 8 csatornás mikrofon el˝oer˝osít˝o Kondenzátor mér˝omikrofon Aktív mélynyomó
1030A
Aktív hangsugárzó
PrecisionM60 Cubase SX3
Laptop
Steinberg
4.4.
Audió szoftver
Megjegyzés Érzékenység= 22mV/Pa Frekvenciaátfogás: 19-85Hz Frekvenciaátfogás: 55-18000Hz Rögzítéshez és kijátszáshoz Rögzítéshez és kijátszáshoz
A mérés kiértékelése
A mérés végeztével rendelkezésünkre állnak a gerjesztésre adott válaszjelek az id˝otartományban. Az impulzusválasz meghatározásához a Fourier-transzformáció nyújt segítséget. Áttérve a frekvenciatartományba, a gerjeszt˝o jel és a válaszjel hányadosaként el˝oálló átviteli függvény inverz Fourier-transzformáltjaként megkapjuk az impulzusválaszt. Azonban a mi konkrét esetünkben a mérés célja a modellezés ellen˝orzése, ezért célszeru˝ lehet a vev˝opontok közötti átvitelt vizsgálni. Ennek az az el˝onye, hogy a gerjeszt˝o jelünket közvetít˝o hangsugárzó nem ideális viselkedésének hatását kiküszöbölhetjük. Ugyanis ez a hatás minden vev˝opontot érint, így ha két vev˝opont közötti átvitelt vizsgálok, ez a probléma nem jelentkezik. Fontos megemlíteni, hogy ez a hatás csak a frekvenciamenetb˝ol származó hibákat javítja, a hangsugárzó geometriájából származó hibákat nem. Geometriai hiba alatt itt azt értjük, hogy a modellez˝o eljárásunkban a gerjesztést pontsebesség forrásként definiáljuk, míg a valóságban használt hangsugárzók ennek nem feltétlenül felelnek meg. Ezt a hatást jelen esetben elhanyagoljuk. A kapott eredményeket így összevetve a modell által számított átviteli fügvvényekkel, leellen˝orizhetjük a modellünk helyességét.
39
4.5.
Eredmények
Annak érdekében, hogy eldönthessük, modellezési eljárásunk mennyire fedi a valóságot. A 4.5, a 4.6 és a 4.7 ábrákat vizsgáljuk. Ezeken az ábrákon kiválasztott mikrofonpárok közti mért átvitelthasonlíthatjuk össze a perturbáció módszer által számított átvitellel. A következ˝o észrevételeket tehetjük. El˝oször a mérési jel eleje és vége tunik ˝ fel, ez a zajos rész kiválóan megmutatja nekünk a gerjesztésként használt hangsugárzó frekvenciaátfogását, ugyanis azokon a tartományokon azért kaptunk ilyen zajos jelet, mert a hangsugárzó képtelen volt megfelel˝o mennyiségu˝ energiát táplálni a rendszerbe. Megfigyelhetjük továbbá, hogy a csúcsok kicsit eltolódnak oldalra, ez annak az eredménye, hogy a szoba méretei nem teljes pontossággal lettek beadva a modellez˝o programba, így a módushelyek kissé eltolódhattak. Látható továbbá, hogy az alacsonyabb frekvenciákon nagyon nagy pontossággal illeszkednek a görbék, azonban magasabb frekvenciákon csupán jellegében követi a modell a mérés eredményeit. Ennek oka abban keresend˝o, hogy míg alacsony frekvencián jól becsültük a falelnyelés mértékét, magasabb frekvenciákon ez a becslés már nem kielégít˝o. Frekvenciafügg˝o elnyelés bevezetésével ez a probléma valószínuleg ˝ nagy hatékonysággal megoldható. 20
10
Amplitúdó[dB]
0
−10
−20
−30
−40 0
20
40
60 80 Frekvencia[Hz]
100
120
4.5. ábra. A 16-os és 57-es vev˝opontok közti átvitel. Piros: perturbáció módszerrel a várható érték, kék: mérési eredmény.
40
15 10
Amplitúdó[dB]
5 0 −5 −10 −15 −20 30
40
50 Frekvencia[Hz]
60
70
4.6. ábra. A 37-es és 57-es vev˝opontok közti átvitel. Piros: perturbáció módszerrel a várható érték, kék: mérési eredmény.
20 15 10
Amplitúdó[dB]
5 0 −5 −10 −15 −20 −25 −30 −35 0
20
40
60 80 Frekvencia[Hz]
100
120
4.7. ábra. A 44-es és 57-es vev˝opontok közti átvitel. Piros: perturbáció módszerrel a várható érték, kék: mérési eredmény.
41
5. fejezet Összefoglalás Dolgozatunk célja egy akusztikai térszámító program megalkotása volt, amely képes egy véletlen paraméter, a véletlen fal elnyelés kezelésére. Dolgozatunkban definiáltuk az akusztikai térszámítás feladatát és bemutattuk annak determinisztikus megoldását, a spektrális végeselem módszert, illetve az azon alapuló akusztikai térszámító szoftvert. Egyszeru˝ példákon keresztül ismertettünk két sztochasztikus rendeszermodellez˝o eljárást, a perturbáció módszert és a Monte Carlo módszert. Ezután megmutattuk, hogy az eredeti determinisztikus térszámító programot milyen átalakításokkal tettük alakalmassá véletlen paraméterek kezelésére. Megalkottuk egy létez˝o terem modelljét, majd a teremben történt hitelesít˝o mérés segítségével megvizsgáltuk eredményeink helyességét. Összességében elmondhatjuk, hogy a mérés alátámasztotta a modellünk és a számítási módszer helyességét. Egyben azt is megmutattuk, hogy amennyiben arra lehet˝oség adódik, sztochasztikus modellez˝o eljárásoknál érdemes keresni alternatívákat a Monte Carlo módszer mellé, jelen esetben a perturbáció módszerét, mert ezzel a számításaink ideje drasztikusan lecsökken, így a modellezési rugalmasság ugyanilyen arányban megn˝o. Sikerült megalkotnunk egy valószínuségi ˝ paraméter, a véletlen falelnyeléssel dolgozó programot, amely hatékony akusztikai térszámítást tesz lehet˝ové. Továbblépési lehet˝oség lehet egy olyan szoftver fejlesztése, amely képes számításokat végezni olyan esetekben, ahol a falak elnyelése nem azonos a fal teljes felületén, illetve képes az elnyelés frekvenciafüggését figyelembe venni. Ezen kívül a forrás pontosabb, nem pontforrásként való modellezése szintén továbbfejlesztési irány.
42
Irodalomjegyzék [1] Ivan T. Dimov. Monte Carlo Methods for Applied Scientists. World Scientific, London, 2008. [2] Stijn Donders, Dirk Vandepitte, J. Van de Peer, and Wim Desmet. Assessment of uncertainty on structural dynamic responses with the short transformation method. Journal of Sound and Vibration, 288:523– 549, 2005. [3] P. Fiala. Puta toolbox for matlab. user’s manual. BME Híradástechnikai Tanszék, 2007. [4] P. Fiala. Development of a numerical model for the prediction of groundborne noise and vibration in buildings. PhD thesis, Budapesti Muszaki ˝ és gazdaságtudományi Egyetem, 2008. [5] Roger G. Ghanem and Pol D. Spanos. Stochastic Finite Elements. A spectral approach. Dover Publications, 2003. [6] Cs. Huszty and F. Augusztinovicz. A teremakusztikai méréstechnikáról és a pásztázó szinuszjelekr˝ol. Akusztikai Szemle, 8(1–2):29–41, 2008. [7] Cs. Huszty and Attila B. Nagy. BMEVIHIAV67 teremakusztika. Hallgatói segédlet, BME Híradástechnikai Tanszék, 2008. [8] Allan D. Pierce. Acoustics. An Introduction to its Physical Principles and Applications. McGraw Hill, New York, 1981. [9] A. Rényi. Valószínuségszámítás. ˝ Tankönyvkiadó, Budapest, 1966. [10] M. Schevenels. The impact of uncertain dynamic soil characteristics on the prediction of ground vibrations. PhD thesis, Catholic University Leuven, 2007.
43
[11] O.C. Zienkiewicz. The finite element method. McGraw-Hill, third edition, 1986.
44