Akoestische modellering van axisymmetrische openingen met Rayleigh elementen

Akoestische modellering van axisymmetrische openingen met

Rayleigh elementen

Citation for published version (APA):

Starreveld, J. P. (1997). Akoestische modellering van axisymmetrische openingen met Rayleigh elementen. (Herziene dr. redactie) (DCT rapporten; Vol. 1997.054). Technische Universiteit Eindhoven.

Document status and date: Gepubliceerd: 01/01/1997

Document Version:

Uitgevers PDF, ook bekend als Version of Record

Please check the document version of this publication:

• A submitted manuscript is the version of the article upon submission and before peer-review. There can be important differences between the submitted version and the official published version of record. People interested in the research are advised to contact the author for the final version of the publication, or visit the DOI to the publisher's website.

• The final author version and the galley proof are versions of the publication after peer review.

• The final published version features the final layout of the paper including the volume, issue and page numbers.

Link to publication

General rights

Copyright and moral rights for the publications made accessible in the public portal are retained by the authors and/or other copyright owners and it is a condition of accessing publications that users recognise and abide by the legal requirements associated with these rights. • Users may download and print one copy of any publication from the public portal for the purpose of private study or research. • You may not further distribute the material or use it for any profit-making activity or commercial gain

• You may freely distribute the URL identifying the publication in the public portal.

If the publication is distributed under the terms of Article 25fa of the Dutch Copyright Act, indicated by the “Taverne” license above, please follow below link for the End User Agreement:

www.tue.nl/taverne

Take down policy

If you believe that this document breaches copyright please contact us at: openaccess@tue.nl

providing details and we will investigate your claim.

Akoestische

modeUeri

axisymmetrische

Ray

le ig h

ele

me

nlten

J.P.

Starreveld 420726

Skage-Rapport 97.054 WFW-dynamica Eindhoven, Juli 1997 (herziene druk)

Begeleid er: ia.. A.H.W.U. KGjpers Afstudeerdocent: prof.dr.ir. D.H. van Campen

Technische

Universiteit Eindhoven

Inhsudsopgave

Samenvatting 1 1 Inleiding 1 1.1 Probleemstelling

. . .

1 1.2 Opdrachtomschrijving

. . .

1 1.3 Gevolgde strategie

. . .

1

2 Theoretische afleiding van de geluidsafstraling van een baffle 3

2.1 Rayleigh integraal

. . .

3 2.2 Rayleigh integraal voor axisymmetrische constructies met niet-axisymmetrische

randvoorwaarden

. . .

5 Constructie van de benaderingsoplossing voor de Rayleigh integraal

. . .

7

2.3.1 7

2.3

Berekening van de Rayleigh integraal over het contactvlak

. . .

2.3.2 Berekenen van de integraal H ,

. . .

8

3 Koppeling van Rayleigh elementen aan BEM 11

3.1 Theoretische afleiding ten behoeve van koppeling

. . .

12 Koppelen van axisymmetrische constructies met niet-axisymmetrische rand- voorwaarden

. . .

12 3.2

4 Numerieke implementatie 15

4.1 Programma opbouw

. . .

15 4.2 Numerieke integratie

. . .

17 4.2.1 Berekenen van de integraal

H ,

. . .

17

Berekening van de Rayleigh integraal over het contactvlak

. . .

4.3 Koppeling van Rayleigh integraal aan BEM

. . .

20

iv

5 Verificatie met analytische oplossingen 23

5.1 Oscillerende zuiger in een vlakke oneindige plaat

. . .

23

Geluidsdruk op het wandoppervlak voor laag frequente oscillaties

. . .

5.1.2 Geluidsdruk in het verre veld voor een zeer locale wandtrilling

. . . .

25

5.1.3 Geluidsdruk op de as van symmetrie

. . .

26

5.1.4 De ontwikkeling van ‘druklobben’ bij toenemende axiale afstand

. . .

27

5.2 Geluidsafstraling van een pijp met een opening aan één zijde

. . .

27

5.3 Pulserende cilinderwand in een akoestisch baffle

. . .

29

5.1.1 23 6 Conclusies en aanbevelingen 33 6.1 Conclusies

. . .

33

6.2 Aanbevelingen

. . .

34

A Opbouw van het BEM-pakket (MATLAB-procedures) 37 B Numerieke integratie 39 B.l Formule van Euler (trapeziumregel)

. . .

39

B.2 Formule van Simpson

. . .

39

B.3 QUAD8 . standaard MATLAB-procedure

. . .

40

B.4 GUASSHM . BARD procedure

. . .

40

C Waarden van [ en wint ten behoeve van Gauss-Legendre-integratie 41 D Log-Gauss-transformatie 43 D.1 Geval 1

. . .

43

D.2 Geval2

. . .

44

D.3 Geval3

. . .

45

Hoofdstuk

1

Inleiding

1

.I

Probleemstelling

Het project ‘Het geluidsbewust ontwerpen van MRI-systemen’ heeft als doelstelling om de mechanica en akoestiek van ‘Magnetic Resonance Imaging’ (MR1)-scanners beter t e begrijpen en t e verbeteren. De in dit project opgedane kennis moet worden gebruikt om al in de ontwerpfase van deze scanner gericht te kunnen zoeken naar een zo stil mogelijk MRI-systeem. In het kader van het bovengenoemde project is in eigen beheer het Boundary Elementen Me- thode (BEM) pakket BARD geschreven, waarmee de geluidsproduktie kan worden berekend van axisymmetrische constructies veroorzaakt door (niet-axisymmetrische) trillingen.

Voor de modellering van de akoestische afstraling van een MRI-systeem wordt in eerste in- stantie uitgegaan van een zogenaamd gebaffelde cilinder model - een model van de afstraling van een cilinder met aan de uiteinden vlakke oneindige platen (baffles). De beschrijving van die baffles is echter nog niet in het programma aanwezig. Een adequate beschrijving van de akoestische baffle aan de uiteinden van een cilinder moet daarom worden ontwikkeld.

1.2

Opdracht omschrijving

Het doel van de stageopdracht is de ontwikkeling, implementatie en het testen van een (nu- merieke) formulering waarmee de afstraling vanuit de uiteinden van een axisymmetrische con- structie in een akoestische baffle kan worden beschreven. Deze formulering moet gekoppeld worden aan de BEM die reeds geïmplementeerd is voor de berekening van de geluidsafstraling in de cilinder.

1.3

Gevolgde strategie

In hoofdstuk 2 zal een theoretische afleiding volgen van de formulering waarmee de geluids- afstraling vanuit het uiteinde van een axisymmetrische constructie in een akoestische baffle kan worden beschreven. Daarbij is gebruik gemaakt van theorie uit [i] en [2]. Deze afleiding

2 Inleiding

zal leiden tot de zogenaamde ‘Rayleigh-integraal’. Vervolgens zal overgegaan worden tot de constructie van de benaderingsoplossing ten behoeve van de numerieke implementatie. [i], [3] en [4] is de hierbij gebruikte theorie.

Hoofdstuk 3 behandelt de koppeling van de Rayleigh elementen aan de BEM. In [3] en [5] Mioïdt een xethode Seschrever, err, twee afzonderlijke REM sommen o p de contactvlakken t e koppelen. Deze koppeiing wordt verkregen door compatibiiteit i e eisen op de contactu!akken. De implementatie van bovengenoemde formulering is uitgevoerd binnen MATLAB. In hoofd- stuk 4 wordt de opbouw van het programma beschreven. Daarnaast zullen enkele bijzonder- heden binnen het programma aan het licht worden gebracht, zoals de benodigde integratie methoden voor het bepalen van de ‘Rayleigh-integraal’. Ook zal worden ingegaan hoe de numerieke koppeling tussen de twee afzonderlijke BEM sommen t o t stand komt.

In hoofdstuk 5 worden enkele analytische oplossingen van de ‘Rayleigh-integraal’ gepresen- teerd. Ook wordt de oplossing van een gekoppeld probleem gepresenteerd. Deze oplossingen worden vergeleken met de hier ontwikkelde methode.

Hoofdstuk

2

Theoretische afleiding van

de

geluidsafstraling van een baffle

In dit hoofdstuk zal de formulering worden afgeleid waarmee de geluidsafstraling vanuit het uiteinde van een axisymmetrische constructie in een akoestisch bafFle kan worden beschreven.

2.1

Rayleigh integraal

De algemene definitie van geluid is: ‘De mechanische trillingen die zich als golven uitbreiden in elastische media in het gebied van de hoorbare frequenties (16 [Hz] - 16 [kHz])’. De

voortplanting van deze harmonische akoestische golven kunnen worden beschreven met de Helmholtz-vergelij king:

V2#

+

k2# =

o,

waarin

#

de scalaire snelheidspotentiaal die evenredig is met de amplitude van de harmonische drukverstoring

(5

= ikpc# met p de dichtheid) en IC het golfgetal welke gedefinieerd is als:

W

k = -

C ’

waarin w en c respectievelijk de golffrequentie en de geluidsnelheid zijn. Met het theorema volgens Kirchhoff en Helmholtz [2], kan dit ook geschreven worden als:

q(v2

+

k 2 ) # - #(O2

+

k 2 ) @ =

v

*

(W#

- #W),

waarin

KP

een functie van de positie is. Door integratie over het volume

V ,

waarbij alle punten buiten het afstralende oppervlak S en binnen een atmosfeer met een straal

R

bevinden, en vervolgens toepassen van het Gauss-theorema ((V2

+

k 2 ) #

= O binnen

V ) ,

levert:

4 Theoretische afleiding van de aeluidsafstraling: van een baffle

waar IR en

ns

respectievelijk de oppervlakte integraal over het atmosferisch oppervlak en de eenheidsnormaal op het oppervlak S zijn. Het minteken voor de eerste term in het rechterlid wordt veroorzaakt door het feit d a t n S naar buiten gericht is. Door Xi? t e definiëren als zijnde de Green’s functie q ( P , Q ) = e-ikR(P7Q)/R(P, Q )

( R

=

IP-QI)

-een functie die de oplossing is voor een puntbron in een 3 dimensionaal oneindig systeem - wordt het linkerlid gelijk aan C ( P ) + ( P ) . Daarnaast volgt! wanneer

XP

voldoet aan de Sommerfeid straiingsconditie (limR,,

[ R

(g

-

ik+)]

= O), met l i m ~ + ,

XP

= O en met l i m ~ + ,

4

= O d a t

IR

gelijk nul is:

Dit vormt de basis voor de BEM. Daarbij is C ( P ) = 47r wanneer P in het vrije veld ligt en

C ( P )

= 2n als P op het oppervlak ligt. De algemene definitie van C ( P ) luidt:

C ( P )

= 4n

+

l’ls

&

(&)

WQ).

De modellering van de gebaffelde cilinder wordt als een axisymmetrisch probleem beschouwd en is schematisch weergegeven in figuur 2.1. Hierbij wordt een cilindrisch coördinaten stelsel

Z

Figuur 2.1: Schematische weergave van de gebafieide cilinder

gehanteerd. Ten gevolge van de in axiale richting volledig reflecterende wand kunnen we het zogenaamde spiegeleffect toepassen. Het spiegelpunt van P is in de figuur aangegeven met

P r .

Het punt Q ligt op het contactvlak tussen cilinder en baffle. Door het spiegeleffect toe t e passen voor het punt Q kan de Green’s functie geschreven worden als:

2.2 Rayleigh integraal voor axisymmetrische constructies met

niet-axisymmetrische randvoorwaarden 5

waarbij :

r p ; r p l ; TQ

z p ; z p r ; z ~ f axiale afstand van respectievelijk P ,

P‘

en Q tot de wand

6 = @Q - 6 p

radiale afstand van respectievelijk P , P’ en Q t o t de rotatie as

hoek tussen Q en

P

Omdat voor een gebaffelde cilinder r p = r p ~ , z p = - Z ~ I en ZQ = O geldt, kun d e Green’s

functie geschreven worden als:

R(P, Q ) = [ r p

+

r& - 2rprQ cos O

+

z;]$.

De factor 2 in

Q(P,

Q) kan fysisch worden verklaard door het wandoppervlak te beschouwen als een infinitesimale dunne permeabele plaat. De normaalsnelheid - is dan aan beide zijde van de plaat van gelijke grootte en richting (het medium links van de plaat beweegt naar links en het medium rechts van de plaat beweegt ook naar links). Vandaar d a t Q(P,Q) voor vermenigvuldigd wordt met de waarde 2 ((2.2)).

Vanwege het feit d a t Q f Q ( z Q ) kan de afgeleide van de Green’s functie geschreven worden als:

82,

Beschouwen we weer de infinitesimale dunne permeabele plaat

,

dan kan deze formule samen met (2.2) worden verklaard. Wanneer ter plaatse P en P‘ een gelijke druk N + ( P ) wordt

opgewekt, dan wordt in het punt Q links en rechts van de plaat een gelijke druk waargenomen. De netto bijdrage van de oppervlakte druk in de Kirchhoff-Helmholtz integraal is daarom gelijk aan nul.

Vandaar d a t de integraal (2.1) nu overgaat in:

Deze vergelijking wordt ook wel de Rayleigh integraal genoemd. Omdat

C ( P )

geen functie van ZQ is, geldt altijd d a t

C ( P )

= 47r. Dit volgt uit de algemene definitie van C ( P ) (pagina 4).

2.2

Rayleigh integraal voor axisymmetrische constructies met

niet-axisymmetrische randvoorwaarden

Het

BEM

pakket BARD is gebaseerd op een theorie die de geluidsafstraling berekend voor axisymmetrische constructies veroorzaakt door niet-axisymmetrische trillingen. Vandaar dat een dergelijke theorie ook moet worden toegepast voor de afstraling vanuit de uiteinden van een cilindische constructie in een akoestische baffle.

, - i W P , & )

Ten gevolge van d e axisymmetrie kunnen de niet-axisymmetrisch variabelen

4(f)

,

R(P,Q)

en -(Q) 84 worden uitgedrukt in een Fourier-reeks over een omwentelingshoek met een periode

6 Theoretische afleiding van de geluidsafstraling van een baffle

van 2n volgens:

00

4 ( P ) =

2

4 0 (PI

+

C(ó”

(

P )

sin 18p

+

@ ( P ) cos ZOp)

1=1 waarin: 1 1

I{:

= - sin m8pHm

Km

= - cos mop H,, ?r Ir met:

Door de introductie van deze Fourier-reeksen kan (2.3) nu worden geschreven als (met dS(Q) =

7.Q d 8 ~ dLQ):

Rangschikken van termen en met gebruikmaking van het feit d a t de integralen over de termen sin m 8 ~ sin n 8 ~

,

sin m 8 ~ cos n8Q en cos m6Q cos n 8 ~ slechts ongelijk nul zijn wanneer m = n,

leidt n a integratie over de hoek OQ tot:

2.3 Constructie van de benaderingsoplossing voor de Rayleigh integraal 7

Uit deze vergelijkingen dienen voor elk gewenst punt P , #o(P), & ( P ) en

&(I')

te worden bepaald, wanneer +O'(&),

4E'(Q)

en &'(Q) bekend zijn, of juist andersom.

Tot dusver is er nog geen sprake van een benadering.

2.3 Constructie van de benacieringsöpiossing voor

de Ray-kigh

integraal

Een eerste benadering van (2.4) is het meenemen van slechts M Fourier-termen in plaats van een oneindig aantal.

In de BEM som zijn we vooral geïnteresseerd in de c,b(P)'s en de "-(P)'s op het contactvlak. Deze zijn nodig voor de koppeling van het gebaffelde cilinder model aan BARD (zie hoofdstuk 3). Daartoe verdelen we het contactvlak op in N elementen. Daarbij wordt gewoonlijk gebruik gemaakt van kwadratische lijn-elementen ([i]). Met dit 3-knoopselement zal (2.4) worden gediscretiseerd.

82,

2.3.1 Berekening van de Rayleigh integraal over het contactvlak

In paragraaf 2.2 hebben we gezien dat de Rayleigh integraal voor een axisymmetrische con- structie met niet-axisymmetrische randvoorwaarden wordt gegeven door (2.4) :

4 r 4 o ( P ) = -2 ~ ~ ( Q ) H o ~ Q ~ L Q 47r4:(~) = -2

s,

s,

s,

4n'(Q)H,rQ dLQ voor n = i , .

.

. )

Ad

4 r 4 i ( p ) = -2 _{&'(&)H,TQ dLQ} voor n = 1,.

.

.

,

M .

Zoals zojuist vermeldt, wordt het contactvlak verdeelt in N elementen. We concentreren ons nu even op één element. Binnen dit element introduceren we het locale coördinatensysteem zoals weergegeven in figuur 2.2. Binnen dit element kan voor respectievelijk straal r ( < ) en

1 2 3

e

I ~

a

r=rl r=r2 C=l r=r3 5 - 1

15

Figuur 2.2: Locale coördinatensysteem Jacobiaan J ,

(<)

geschreven worden:

8 Theoretische afleiding van de geluidsafstraling van een baffle

Op dezelfde wijze kunnen nu ook de overige variabelen uit (2.5) worden uitgeschreven (ele-

Door nu de gedefinieerde kolommen te substitueren in (2.5) levert dit na sommatie over alle elementen:

J

d t ] voor n = 1,.

.

.,

M (2.6)

voor n = 1,.

.

. , M ,

met:

HieïUq dienen de integraler, ~ u m e r i e k te werden bepaald. Zoals nader in hoefdstiik 4 zal worden besproken, dient hier een Gauss-Legendre- en Log-Gauss-integratie schema t e worden gehanteerd.

2.3.2 Berekenen van de integraal

H,

Zoals in paragraaf 2.2 werd beschreven, wordt

H ,

gegeven door:

2.3 Constructie van de benaderingsoplossing voor de Rayleigh integraal 9

Hierin kunnen weer de gedefinieerde kolommen uit d e vorige paragraaf worden gesubstitueerd, waardoor deze gebruikt kan worden in combinatie met (2.6).

Bij het oplossen van deze integraal kan bijvoorbeeld gebruik worden gemaakt van de formule van Euler (trapezium regel), de formule van Simpson of een hogere orde integratie schema. Over welke methode het beste kaiaI; worden toegepast, w d t ir, hoofdstuk 4 nader ingegaan.

Hoofdstuk

3

Koppeling

van Rayleigh elementen

aan

BEM

M e t behulp van standaard BEM, (2.1):

C ( P M P ) =

/’

( ~ V W P ,

Q )

-

w,

Q)WQ))

.nS

WQ),

S

is het in principe mogelijk o m de akoestische afstraling van een MRI-systeem te modelleren. Tenminste in het geval dat het golfgetal k in domein 1 en 2 gelijk zijn (zie figuur 3.1). Het op

‘t

2

. P

0

‘11s

. P

o

N elementen

Figuur 3.1: Wijze van modellering van het MRI-systeem

deze manier modelleren heeft echter twee nadelen:

1. In elk punt

P

dient bovengenoemde integraal geïntegreerd te worden over het oppervlak

S

i

en 5’2.

2. Men kan zich afvragen, hoe men de vlakke oneindige plaat kan modelleren.

Dit hoofdstuk zal de koppeling van de Rayleigh elementen aan het BEM pakket BARD behandelen, waardoor de zojuist genoemde nadelen komen te vervallen.

12 Koppeling van Rayleigh elementen aan BEM

3.1

Theoretische afleiding ten behoeve van koppeling

In plaats van slechts één integraalvergelijking uit t e gaan kan het akoestisch probleem ook worden opgedeeld in tweeën. Eén voor de modellering van de geluidsafstraling binnen de cilinder en één voor de modellering van de geluidsafstraling naar de vrije ruimte.

De integraalvergelijking voor de akoestische afstraling binnen de cilinder - domein 1 - wordt

beschreven middels:

Deze geluidsproduktie kan worden berekend met behulp van het BEM pakket BARD. Voor de geluidsafstraling naar de vrije ruimte wordt uitgegaan van de in hoofdstuk 2 afgeleide Rayleigh integraal (2.3):

O p SI gelden de voorgeschreven randvoorwaarden, terwijl op s compatibiliteit wordt geëist. Dat wil zeggen dat op het contantvlak s voor de continuïteit van de normaalsnelheden wordt geëist dat:

en voor de continuïteit van de druk wordt geëist dat:

p i + i ( Q ) = p242(Q) OP S.

3.2 Koppelen van axisymmetrische constructies met niet-axi-

symmetrische randvoorwaarden

De binnen het BEM pakket BARD aanwezige integraalvergelijkingen zijn zoals gezegd geba- seerd op een theorie die de geluidsafstraling berekend voor axisymmetrische constructies ver- oorzaakt door niet-axisymmetrische trillingen. Ter beschrijving van deze niet-axisymmetrische trillingen zijn net zoals in hoofdstuk 2 Fourier-reeksen ingevoerd. Daarbij worden weer slechts

M

Fourier-coëfficiënten meegenomen. Ook wordt de lijn LI en 1 binnen domein 1 opgedeeld in NI

+

N kwadratische lijn-elementen. Voor domein 1 komen de integraalvergelijkingen er dan als volgt uit de te zien ([i]):

C(p)4Ol(p)

Ni +N 1

[ 6 û m i T ~ ~ N ( t ) H ó ( t ) T Q ( t ) J m ( t ) d t - $ûmllT

LI

N (t)

H O ( e ) ‘Q

(t)

Jm

(e)

d t ]

m=l

3.2 Koppelen van axisymmetrische constructies met niet-axisymmetrische

randvoorwaarden 13

met:

en voor de Rayleigh integraal gelden de vergelijkingen (2.6). Deze vergelijkingen samen met de compatibiliteitseisen maken vervolgens de koppeling mogelijk. Op de numerieke imple- mentatie zal in hoofdstuk 4 uitvoerig worden ingegaan.

Hoofdstuk

4

Numerieke implementatie

In dit hoofdstuk zal besproken worden hoe het boundary elementen methode programma is opge- bouwd. Daarnaast wordt ingegaan op de benodigde numerieke integratie-methoden. Ook wordt

de numerieke koppeling van beide BEM sommen beschreven.

4.1

Programma opbouw

Voor de manier van opbouw is gekeken naar a n c x e programma’s. Neem bijvoorbec,,, de programmatuur beschreven in [ 3 ] , de .mud file van het EEM pakket MARC, de invoerfile van ANSYS en de programmatuur behorende bij het vak: ’Applied Computational Mechanics’. Het hier beschreven BEM pakket is globaal als volgt opgebouwd:

1. Het opvragen of creëren van een mesh met randvoorwaarden en deze weergeven op het scherm. De mogelijke randvoorwaarden zijn: $0, $;,

4:,

4:’’

4;

en/of

4:

voorschrijven voor n Fourier-coëfficiënten.

2. Het invoeren van het golfgetal k en de materiaaleigenschappen (binnen het hoofdpro- gramma (RAYLE1GH.M) kan het aantal integratie punten ten behoeve van Gauss- Legendre- en Log-Gauss-integratie worden opgegeven en kan de tolerantie van de inte- gratie van H , worden opgegeven].

3 . Het bepalen van de Rayleigh integralen voor de elementen en deze assembleren tot de

zogenaamde matrices

G

en H - zie paragraaf 4.3 (de integralen in de knooppunten Q).

4. Oplossen van het stelsel vergelijkingen:

H(P,, Q ,

&)$(&I

= G(P,, Q ,

Q ) ~ Y Q )

of het koppel probleem Dit levert de onbekende $0, q!&,

+;,

&’, 4;

en/of

4:

op de elementen op.

5. Het bepalen van de Rayleigh integralen voor de veldpunten en deze assembleren tot de matrices

G

en H (de integralen in de veldpunten

P ) .

6. Oplossen van het stelsel vergelijkingen:

16 Numerieke implementatie

Dit levert de onbekende 4 0 ,

4;,

4;,

&’,

4;

en/of

4;

in de veldpunten op.

In figuur 4.1 staat de opbouw nog wat meer in detail weergegeven. Hierin is duidelijk het

:

I

Plot axisymmetrische mesh

I

Ko

k.

el

probleem

Plot randvoo waarde n

+

Invoeren van materiaal eigenschappen

Bepalen van Rayleigh-integralen voor rand-elementen

&

Assembleren tot matrices

1

1

Bepalen van drukken en snelheden op de rand-

I

Vraag een bArd-mesh op

I

Plot contact-elementen

I

voor contact-elementen

1

Bepalen van drukken en snelheden op de contact- elementen en op de buiswand

voor veld-punten

Bepalen van drukken en snelheden in de veld-punten

Figuur 4.1: Opbouw van het BEM pakket

4.2 Numerieke integratie 17 0.6 0.5 0.4 - _c I ä 0.3- o o I 0.2- L 0.1 O -0.1

bijlage A staat hetzelfde schema nogmaals gegeven, maar nu is hierin aangegeven welke pro- cedures achtereenvolgens worden aangeroepen.

- - - Simpcon n=24 iniemailen Quad8 tolerantie = i e 3 - - I I / ) I (I I I - ' I 1 1 ' I I / I . I . . --.I

.--

' ' ' ' ' ' in0

4.2

Numerieke integratie

In paragraaf 2.2 is reeds ter sprake gekomen, d a t voor het oplossen van de Rayleigh ifite- graal diverse integratie-methoden moeten worden toegepast. In deze paragraaf zullen deze methoden worden besproken.

4.2.1 Berekenen van de integraal

H ,

Bij het oplossen van de integraal:

kan gebruik worden gemaakt van de formule van Euler (trapezium regel), de formule van Simpson of een hogere orde integratie methode, zoals QUAD8 - een standaard MATLAB-

procedure - en GAUSSHM - een procedure die onder andere gehanteerd wordt binnen BARD.

Deze integratie methoden worden in bijlage B kort besproken.

In figuur 4.2 zijn de resultaten van de absolute waarde van H , ( n = O) weergegeven als functie van T-Q ( r p = i). In de linker grafiek staan de resultaten weergegeven wanneer voor de Euler- en Simpson-integratie 24 intervallen worden genomen. In de rechter grafiek staat de relatieve fout van (H,( voor de diverse methoden ten opzichte van GAUSSHM weergegeven. Zoals uit

1 00

rQ

rQ

Figuur 4.2: Integratie van H , met Euler, Simpson, QUAD8 en GAUSSHM

de grafieken blijkt, treedt bij T-Q = 1 een singulariteit o p ( R ( P , Q ) = IQ -

PI

= O). Hierbij

18 Numerieke implementatie

Daarnaast is waar te nemen d a t de resultaten met behulp van de formule van Euler en Simpson afwijken van de resultaten met de hogere orde methoden. Een toename van het aantal intervallen blijkt echter wel te leiden tot een afname van de relatieve fout. Doordat de procedure QUAD8 de mogelijkheid biedt tot het opgeven van een tolerantie, is deze bij de verdere numerieke implementatie toegepast.

4.2.2 Berekening van de Rayleigh integraal over het contactvlak

Zoals uit de vorige paragraaf blijkt treedt er een k-singulatiteit op. De singulariteit treedt alleen op in het element met de node waarvoor

R(P,

Q ) = IQ - PI = O. Met deze wetenschap kunnen we nu twee gevallen onderscheiden:

1. integratie over elementen zonder &-singulariteit 2. integratie over elementen met &singulariteit

Volgens [3] is het daarbij gebruikelijk om voor geval 1 Gauss-Legendre-integratie toe t e pas- sen. Voor geval 2 dient echter een nauwkeurigere integratie methode te worden toegepast, zoals een hogere orde integratie, of een speciale formule (bijvoorbeeld Log-Gauss). Hier is gekozen voor Log-Gauss, omdat het reële deel van Hn zich in de buurt van de singulariteit logaritmisch gedraagt.

Opgemerkt dient t e worden dat in plaats van Log-Gauss-integratie ook gewoon Gauss-Legendre- integratie kan worden toegepast. Echter het grote voordeel van Log-Gauss is d a t de rekentijd sterk wordt gereduceerd bij gelijk blijvende nauwkeurigheid (minder aantal integratie-punten nodig).

Gauss-Legendre-int egrat ie

Gauss-Legendre-integratie is niets anders dan het invullen van een waarde (hier

E )

in de t e integreren vergelijking, gevolgt door een vermenigvuldiging met een weegfactor wint. De waarde vac

E

e n wint zijr, af'hankelijk van het aantal integratie-punten. De Rayleigh integraal (2.6) komt er dan als volgt uit t e zien:

m=i

L

ip=l

Deze integratie is terug t e vinden in de procedure GEXT.M.

In bijlage C is een tabel met waarden van

E

en wint gegeven voor een aantal integratie punten van 1 tot en met 12. Dit staat ook in de procedure GAUSS.M.

4.2 Numerieke integratie 19

Log- Gauss-integratie

De t e integreren termen in de Rayleigh integraal luiden voor een kwadratisch element:

Hierbij treedt de singulariteit op in

H,,

zoals we al eerder in deze paragraaf zagen. Wanneer

P

Q dicht nadert kan het reële deel van Hn beschouwd worden als een functie die evenredig is met In

(:)

([3]). De integratie grenzen van deze functie lopen van O tot 1, in plaats van -1 tot 1.

Log-Gauss-integratie is van de vorm:

Deze t e integreren functie lijkt sterk op

si([).

Vandaar d a t we deze Log-Gauss-integratie methode willen transformeren tot in de vorm van de Rayleigh integraal. Daarbij zijn de in figuur 4.3 aangegeven gevallen t e onderscheiden. In bijlage D zijn de transformaties voor deze drie gevallen uitgevoerd. In tabel 4.1 staan de

geval 1 geval 2 geval 3

I

-

I

5 - 1 5=0 $=i 5'+0 5=1 5=-1 g=0 k=l

q=O q-1 q =1

?i

=q=o q=l q=l q = O

.=knooppunt X=singulariteit

Figuur 4.3: Drie gevallen van Log-Gauss-integratie

resultaten. Wanneer nu de functies

si([)

voor de drie verschillende gevallen grafisch worden

Tabel 4.1: Log-Gauss-transformaties

weergeven, zien we d a t - ten gevolge van het gedrag van de vormfuncties - de Log-Gauss-

integratie alleen mag worden toegepast op het reële deel van

Hn

in het knooppunt waarin de singulariteit optreedt (zie figuur 4.4). Voor de overige knooppunten en het imaginaire deel van H, dient gewoon Gauss-Legendre-integratie t e worden toegepast. Deze integraties zijn terug

20 Numerieke implementatie Singulariteit in knooppunt 1 61 I Singulariteit in knwppunt 2 xi -2

t

-16' I -1 -0.8 -0.6 -0.4 -0.2 O 0.2 0.4 0.6 0.8 1 Singulariteit in knooppunt 3 xi U -1 o -1 -0.8 -0.6 -0.4 -0.2 O 0.2 0.4 0.6 0.8 xi

Figuur 4.4: Grafische weergave van g; voor de drie gevallen

te vinden in de procedure GL0C.M. Standaard worden daarbij 4 integratie punten gebruikt. Het aantal integratie punten voor Gauss-Legendre-integratie kiezen we hieraan gelijk.

In bijlage E zijn de waarden van q en w;,t in een tabei weergegeven voor een aantai integratie punten van 1 to t en met 10. Dit staat ook in de procedure L0GGAUSS.M.

4.3

Koppeling van Rayleigh integraal aan

BEM

uit (2.6) en (3.1) worden berekend. De punten P worden nu zo gekozen dat deze samenvallen met de knooppunten Q uit de domeinen 1 en 2. Op deze manier kunnen de bovenstaande integralen voor elk element in elk punt

P

worden bepaald. We definiëren nu de vectoren

40,

4.3 Koppeling van Rayleigh integraal aan BEM 21

( p p p p ,

,,

n , en

e.

Hierin staan respectievelijk de bekende en onbekende normaalsnelheden

en drukken in de knooppunten

&.

Door de definitie van deze vectoren kunnen de integralen

(9) in de punten P uit domein 1 worden opgeslagen in de matrix G,, enerzijds en de integralen

{ g } in de punten

P

uit domein 2 in de matrix G,, anderzijds. In de matrices H,, en

Hn2

worden de waarden van

C ( P )

en { h } opgeslagen. Op deze manier kan het stelsel vergelijkingen voor domein

1

en 2 geschreven worden als:

Ho, $01 = Go,

01

Ho2& = Go2$& Hn1&1 = Gn14Z1 Hn14;1= ~ n , & 1 Hnz t i 2 = Gn2

{

~ n4 k 2 2 = ~ ne22

{

domein 1 domein 2

Na zorgvuldig rangschikken kunnen de matrices G,, en H,, worden opgesplitst in een gedeelte d a t deel uitmaakt van het oppervlak L I en een gedeelte d a t deel uitmaakt van het contactvlak

1. Met behulp van de in hoofdstuk 3 geïntroduceerde compatibiliteitseisen kan voor iedere

Fourier-coëfficiënt n het volgende stelsel vergelijkingen worden gevonden die moeten worden opgelost:

I

-

-Met behulp van (2.6) respectievelijk (3.1) kunnen met behulp van de nu bepaalde normaal- snelheden en drukken de akoestische afstraling binnen de cilinder en naar de vrije ruimte worden bepaald.

Tot nu toe hebben we ons geconcenteerd op één opening. Het MRI-systeem heeft echter zowel links als rechts een opening. O p dezelfde manier kan hiervoor worden afgeleid:

22 Numerieke implementatie

Deze theorie staat ook binnen de MATLAB-procedure C0UPLE.M. Binnen C 0 U P L E . M worden G,,, G & en Gn, gegeven door B / ( - i k p c ) en worden

H,,, H&

en gegeven door

A . De matrices A en B zijn berekend middels het BEM pakket BARD. Het verschil tussen

de matrix G en B wordt veroorzaakt door een verschil in definitie van de basisvergelijking voor de BEM binnen BARD en RAYLEIGH.

Hoofdstuk

5

Verificatie

met analytische

oplossingen

In dit hoofdstuk zullen enkele analytische oplossingen voor de Rayleigh integraal worden vergeleken met oplossingen verkregen met de hier ontwikkelde BEM. Tevens worden twee voorbeelden gegeven van een gekoppeld probleem.

5.1

Oscillerende zuiger in een vlakke oneindige plaat

Een oscillerende zuiger in een vlakke oneindige plaat is een probleem d a t uitstekend door de Rayleigh integraal beschreven kan worden. Voor een dergelijk probleem bestaan, door het doen van enige aannamen, analytische oplossingen. Deze analytische oplossingen worden gebruikt om het hier ontwikkelde BEM pakket t e controleren.

5.1.1

Geluidsdruk op het wandoppervlak voor laag frequente oscillaties In een vlakke oneindige plaat bevindt zich een zuiger met een straal a van 1 [m]. Deze zuiger osci!!eert met eer, zeer !age f r e y ~ e n t i e (La = 9.01

<

1

[-I).

Met de aamarnen d a t ka

<

1

[-I,

d a t ijn (snelheidsamplitude) constant is over het zuigeroppervlak (Fourier-coëfficiënt O is

gelijk aan 2) en d a t we slechts geïnteresseerd zijn in de geluidsdruk @ o p het wandoppervlak

( 2 = O), kan de Rayleigh integraal worden vereenvoudigd tot ([2], pagina's 218-219):

7T

p =

E.,,

n-

[-ik&

(f)

+

S ( k a ) 2 ]

,

waarin:

De integralen E ( m ) en K ( m ) worden respectievelijk de elliptische integralen van de eerste en tweede soort genoemd en worden gegeven door:

24 Verificatie met analytische oplossingen

Log- Gauss

In figuur 5.1 staat links de mesh weergegeven. Deze mesh bestaat uit 5 3-knoopselementen. Rechts staat de dimensieloze drukamplitude - als functie van voor zowel de analytische

11 4 le-3 1.5 165e- 1 11 2 le-3 8.2722e-1 11 6 le-3 6.0766e-2 11 4 le-6

I.

5070e- 1 6 4 le-3 1.3174e-1 r ! 11- x 10" BEG oplossing a = l [m]

_

Figuur 5.1: Mesh en resultaten van het drukverloop o p het wandoppervlak voor een laag frequent oscillerende zuiger in een vlakke oneindige plaat (La = 0.01

[-I)

als de BEM oplossing gegeven. De BEM oplossing is verkregen met de standaard instellin- gen van RAYLE1GH.M (aantal integratie = 4; tolerantie = le-3). Zoals uit de resultaten blijkt, komt de BEM oplossing uitstekend overeen met de analytische oplossing. Het is echter belangrijk d a t de frequentie f = voldoende klein wordt genomen, omdat de analytische oplossing alleen voor lage frequenties toegepast mag worden.

Uit de resultaten blijkt d a t bij een laag frequente oscillatie van een zuiger in een baffle de geluidsdruk op de as van symmetrie het grootste is. In radiale richting neemt de druk aan de wand af. De verklaring hiervoor is d a t de lucht voor d e zuiger o p de symmetrie-as zich uitsluitend in axiale richting verplaatst, terwijl bij toename van r de lucht steeds meer in radiale richting gaat wegvloeien.

In tabel 5.1 s t a a t de invloed van meshversjning, aar.tal integratie-punten en tolerantie op het verkregen resultaat weergegeven voor zowel Gauss-Legendre- als Log-Gauss-integratie. Zoals

2r

Gauss- Legendre

nodes nint tolerantie relatieve fout

[%I

11 4 le-3 3.1877e-1

11 2 le-3 1.0570e+O

11 6 le-3 1.5449e- 1

11 4 le-6 3.1972e-1

5.1 Oscillerende zuiger in een vlakke oneindige plaat 25

uit de resultaten blijkt, kan, zoals eerder reeds opgemerkt, in de singuliere knooppunten op de elementen het beste Log-Gauss-integratie worden toegepast. Een toename van het aantal integratie-punten leidt voor beide integratie-methoden t o t een verbetering van het resultaat. Meshverfijning leidt niet tot een noemenswaardige verbetering.

5.1.2

Geluidsdruk

in

het verre veld voor een zeer iocaie

wandirilhg

We beschouwen een vlakke oneindige wand waarin zich een oscillerende zuiger bevindt met een straal a van 1 [m]. Deze zuiger trilt met een constante frequentie (ka = 8

[-I)

en heeft een snelheidsamplitude

6,

= 1 [m/s] (Fourier-coëfficiënt O is gelijk aan 2). Wanneer we aannemen d a t de afstand van de zuiger t o t de veldpunten vele malen groter is dan a of

ka2 (bijv. tangentiale afstand = 800a) en de snelheidsamplitude

2,

constant is over het zuigeroppervlak, kan de Rayleigh integraal worden benaderd door

(['LI,

pagina's 225-227) :

p

N f ( 0 , q5)?--1P,

met:

pcûnka2 2Ji (ka sin 0)

ka sin

O

f ( b 4

= -i

J1 wordt wel d e Bessel functie van de eerste orde genoemd en wordt gegeven door:

7r

J I ( q ) =

J,

cos(q cos

4)

sin2

4

d4.

De zuiger wordt gemodelleerd met 5 kwadratische elementen en in het verre veld liggen 30 veldpunten (zie figuur 5.2). In de rechter figuur staat polair met een getrokken lijn de ana-

[ a = l [m] m m x m m f m m m m m m - . m m m f 'M- -.z

Figuur 5.2: Mesh en resultaten van het drukverloop o p grootte tangentiale afstand van een oscillerende zuiger in een vlakke oneindige plaat ( k a = 8

[-I,

tangentiale afstand= 800a) lytische oplossing van de dimesieloze drukamplitude weergegeven en met o de oplossing met behulp van d e

BEM

aangegeven. De resultaten komen zeer goed overeen. Wordt de tangenti- ale afstand echter niet ver genoeg gekozen, dan treden in de buurt van z = O afwijkingen op.

26 Verificatie met analytische oplossingen

Deze afwijkingen worden veroorzaakt door de aannamen die binnen de analytische oplossing zijn gedaan.

Het blijkt dat bij toename van k a het aantal ‘druklobben’ toeneemt. Binnen een ‘lob’ plant de geluidsdruk zich gemakkelijk voort, terwijl daar buiten zich geen geluidsdruk voortplant.

5.1.3

Geluidcdruk

op de as van symmetrie

Beschouwd wordt een vlakke oneindige wand waarin zich een oscillerende zuiger bevindt met een straal a van 1 [m]. Deze zuiger trilt met een constante frequentie ( k a = 1 1 ~

[-I)

en heeft een snelheidsamplitude 8, = 1 [m/s] (Fourier-coëfficiënt O is gelijk aan 2 ) . Wanneer we aannemen dat de veldpunten zich uitsluitend o p de symmetrie-as bevinden, kan de Rayleigh integraal worden vereenvoudigd t o t ( [ 2 ] , pagina’s 232-233) :

Links in figuur 5.3 staat de gebruikte mesh gegeven. In de rechter figuur staat de dimensieloze

Figuur 5.3: Mesh en resultaten van het drukverloop op de symmetrie-as voor een oscillerende zuiger in een vlakke oneindige plaat (ka = 1 1 ~

[-I)

drukamplitude

1

als functie van de axiale afstand f voor de analytische en de BEM oplossing weergegeven. De resultaten komen goed overeen. Alleen dicht bij het oppervlak zien we d a t er een kleine fout (relatieve fout is h2.43

[%I)

optreedt. Door gebruik van meer integratie-punten (standaard 4), kan deze fout worden verkleind. Zo levert bijvoorbeeld 5 integratie-punten een relatieve fout van h0.085

[%l.

(Let op: Bij gebruik van kwadratische elementen kunnen alleen een even aantal integratie-punten worden gekozen.)

Zoals waarneembaar is, blijkt dicht bij het zuiger oppervlak een sterke druk-interferentie op t e treden. In het verre veld echter neemt de axiale druk monotoon af (asymptotisch volgens:

PCl.Unl

5.2 Geluidsafstraling van een pijp met een opening aan één zijde 27

5.1.4

De

ontwikkeling

van

‘druklobben’ bij toenemende

axiale

afstand

Er wordt uitgegaan van een vlakke oneindige wand waarin zich een oscillerende zuiger bevindt met een straal a van 1 [m]. Deze zuiger trilt met een constante frequentie van k a = 20

[-I

en heeft een snelheidsamplitude ijn = 1 [m/s] (Fourier-coëfficiënt O is gelijk aan 2 ) . We zijn nu genteresreed in d e ontwikkeling van de ‘druklobben’ bij toenemende axiale afstand. Deze kan beschreven worden met de volgende analytische oplossing

([a]

pagina’s 236-237) :

met: 1 voor r

<

a H = { O voor r

>

a

R,

= [ ( a - w12

+

z2]f

R~

= [ ( a

+

w)2

+

z213

x =

( & ) ; ( u - W ) .

H wordt de Heaviside- of eenheidstapfunctie genoemd en R, en Ri zijn respectievelijk de kleinste en de grootste afstand van een veldpunt tot de omtrek van de zuiger. AD is de Diffusie integraal en kan met de Fresnel integralen:

worden geschreven als:

AD = - 1 - i e - à n / 2 P {sign(X) - (i - i ) [ C ( X )

+

i S ( X ) ] ) .

2

Links in figuur 5.4 staat de mesh met veldpunten op een axiale afstand van z / a =0.5, 1.0,

1.5 en 2.0 weergegeven. In de rechter figuur staat de dimensieloze drukamplitude

(pLf22,!:,

als functie van de radiale afstand voor de analytische en de BEM oplossing weergegeven. De resultaten zijn uitstekend.

De resultaten geven duidelijk de ontwikkeling van de ‘druklobben’ weer. Een toename van de axiale afstand gaat gepaard met een verschuiving van de ‘lobben’ in radiale richting.

5.2

Geluidsafstraling van een pijp met een opening aan één

zijde

Aan de ene zijde van een pijp bevindt zich een oscillerende zuiger en aan de andere zijde een opening (zie figuur 5.5). De binnenstraal van de pijp is gelijk aan de straal van de zuiger en

28 Verificatie met analytische oplossingen I m Y m m m m m m m m m m m m m m m m x m m m m m m m m m m m m Y m m m m m m m m a m m m m m m m x a = b i m m m m m m x x m m m m m x m m m m m m m m m m m

-

rla

Figuur 5.4: Mesh en resultaten van de visualisatie van de ontwikkeling van ‘druklobben’ bij toenemende axiale afstand voor een oscillerende zuiger in een vlakke oneindige plaat ( k a = 20

[-I

)

Figuur 5.5: Pijp met aan de ene zijde een oscillerende zuiger en aan de andere zijde een opening

heeft een straal a van 1 [m]. De buitenstraal van de pijp is 2a en de pijplengte

I

is 3a. De zuiger

oscilleert met een constante dimensieloze frequentie ku van 0.85

[-l.

De snelheidsamplitude

ijUn is 1 [m/s] (Fourier-coëfficiënt O is gelijk aan 2). In het artikel [5] zijn resultaten voor een dergelijk probleem gegeven. Deze zijn verkregen op basis van een BEM d a t het volledige probleem (‘enkel-domein’) kan modelleren.

Met het BEM pakket BARD en de hier ontwikkelde BEM is een koppeling voor dit probleem uitgevoerd. De pijp is gemodelleerd binnen BARD wat de matrices A en B en een mesh ople- vert. Aan de ene zijde van de pijp wordt de constante snelheidamplitude ijn voorgeschreven, terwijl de andere zijde wordt gekoppeld aan een baffle. Het beschouwen van het uiteinde van de pijp als zijnde een baffle is een benadering.

5.3 Pulserende cilinderwand in een akoestisch baffle 29

90 - gekoppelde EEM oplossing (banle)

% oplossing cilindrische pijp

60

- gekoppelde EEM oplossing (banle)

% oplossing cilindrische pijp

x

x o

zia

Figuur 5.6: Resultaten van de geluidsdruk binnen de pijp en de afstraling naar de vrije ruimte

( k a = 0.85

[-I;

tangentiale afstand = l o a )

drukamplitude

&

op de pijpwand als functie van de plaats weergegeven. Rechts is polair de drukamplitude op een tangentiale afstand van 10a weergegeven. Er treden behoorlijke verschillen o p in de resultaten. Door het gebruik van ‘Rayleigh-elementen’ treedt er meer reflectie op, dan bij de pijp-opening met een buitenstraal van 2u. Hierdoor is de druk aan de buitenwand voor het gekoppelde probleem hoger. Deze extra reflectie leidt tevens tot een verandering van de drukamplitude binnen de pijp. Dit verklaart de verschillen in de resultaten.

5.3

Pulserende cilinderwand in een akoestisch baffle

Binnen BARD is een optie (‘BaffledDuct’) aanwezig om de geluidsafstraling binnen een ge- baffeide cilinder t e modelleren met niet-axisymmetrische randvoorwaarden. Met behulp van deze optie zal het gebruik van niet-axisymmetrische randvoorwaarden binnen de koppeling tilsseri. het, REM pakket, RARD en de hier ontwikkelde REM wsrden gecmtrdeerd.

Daartoe wordt binnen BARD de cilinder gemodelleerd. Dit levert de vereiste matrices A en B en een mesh. Vervolgens worden beide zijden gekoppeld aan de ‘Rayleigh-elementen’ en wordt op de cilinderwand de gewenste mode - Fourier-coëfficiënt - voorgeschreven.

Deze procedure is uitgevoerd voor mode O, 1 en 3 voor de frequenties ka = 0.85

[-I

en

k a = 10

[-l.

Voor mode 1 en 3 worden alleen de cosinus-termen van de Fourier-coëfficiënt voorgeschreven. De cilinder heeft een lengte van 3 [m] en heeft een binnenstraal a van 1 [m]. De dichtheid p en de geluidsnelheid c zijn voor de eenvoud gelijk aan 1 gekozen.

In figuur 5.7 staan de resultaten voor de verschillende modes en frequenties weergegeven. Per rij staan respectievelijk de modes 0, 1 en 3 gegeven, terwijl in de afzonderlijke kolommen respectievelijk de frequenties k a = 0.85

[-I

en ku = 10

[-I

staan. Mode 1 is alleen voor een frequentie van ka = 0.85

[-I

uitgevoerd. De resultaten komen voor lage frequenties zeer goed overeen. Voor hogere frequenties zien we d a t bij een toenemend mode-nummer

30 Verificatie met analytische oplossingen O 0 2 - rnode=O; k=0.85 - P (abs) - - Y (abs)

p & v (abs) gek.BEM

rnode=O; k=10 O -1.5 -1 -0.5 O 0.5 1 z -1.5 -1 -0.5 O 0.5 1 1.5 z rnode=3: k=10

p & v (abc) gek.BEM

0.81

-1.5 -1 -0.5 O 0.5 1

2

5

Figuur 5.7: Resultaten van d e koppeling tussen BARD en de ‘Rayleigh-elementen’ voor mode O, 1 en 3 en d e frequenties k a = 0.85

[-I

en k a = 10

[-l.

5.3 Pulserende cilinderwand in een akoestisch baffle 31

afwijkingen gaan optreden. Het verhogen van het aantal integratie-punten en/of de tolerantie leidt daarbij niet tot noemenswaardige verbeteringen. We zien ook d a t bij toename van

k en n de periodetijd van een trilling afneemt, waardoor ook de golflengte afneemt. Een afname van de golflengte bij constant blijvende elementlengte leidt tot het minder nauwkeurig kunnen beschrijven van een trilling. Zo geldt er de regel: ‘Er zijn minstens 6 knooppunten per golflengte noodzakelijk om een trilling goed te kunnen beschrijven’. Meshverfijning leidt echter alleen t o t een visuele verbetering en niet tot het kleiner worden van de afwijkingen. Bij de resultaten die verkregen zijn met de BARD-optie ‘BaffledDuct’, is een trend waarneembaar d a t bij toename van het aantal Gauss-integratie-punten de afwijking kleiner wordt.

Het grote voordeel van de gekoppelde BEM is d a t het de mogelijkheid biedt om zowel de geluidsafstraling binnen de cilinder als naar de vrije ruimte t e beschrijven.

Hoofdstuk

6

Conclusies en aanbevelingen

6.1

Conclusies

Voor de modellering van het akoestisch baffle is uitgegaan van de zogenaamde ‘Rayleigh- integraal’. Voor numerieke implementatie van deze integraal is over gegaan op een BEM- formulering die geschikt is voor axisymmetrische constructies met niet-axisymmetrische rand- voorwaerden. De koppeling is verkregen door compatibiliteit t e eisen op de contactvlakken van beide BEM sommen. Voor de ontwikkelde BEM formulering kan het volgende geconclu- deerd worden:

o Aan de hand van analytische oplossingen is gevalideerd d a t de hier ontwikkelde Boundary

Elementen Methoden (BEM) pakket binnen MATLAB geschikt is om de geluidsafstra- ling vanuit de uiteinden van een cilindrische constructie in een akoestische baffle te beschrijven.

o De ontwikkelde formulering is t e koppelen aan het BEM pakket BARD. Hiermee beschik- ken we over een gereedschap waarmee de geluidsproduktie kan worden berekend voor zowe! binnen als buiten een gebaffelde cilinder veroorzaakt door niet-axisymrnetrische trillingen. Dit model kan mogelijk uitkomst bieden voor het beschrijven van de akoes- tische afstraling van een MRI-systeem.

Aan de hand van de resultaten kan het volgende worden geconcludeerd:

o Het is gebleken dat voor de bepaling van de ‘Rayleigh-integralen’ in de singuliere knoop- punten o p de elementen het beste Log-Gauss-integratie kan worden toegepast.

o Daarnaast blijkt d a t bij een toename van het aantal integratie-punten de nauwkeurigheid van het eindresultaat wordt verbeterd. Ook het verkleinen van de tolerantie leidt tot een verbetering. Daar staat echter tegenover d a t de rekentijd toeneemt. De standaard instelling, n;,t = 4 en tol = l e -

3,

is echter een goed compromis.

o Meshverfijning is zinvol wanneer er minder dan 6 knooppunten per golflengte zijn. Er

34 Conclusies en aanbevelingen

6.2

Aanb evellngen

Het hier ontwikkelde BEM pakket op basis van ‘Rayleigh-elementen’ is in principe een ge- reedschap, d a t zonder enige beperking met BARD kan worden gekoppeld. Toch t o t slot nog enkele aanbevelingen:

o Tijdens d e ontwikkeling van het BEM pakket RAYLE1GH.M is gebleken d a t MATLAB, vanwege d e lange rekentijden, niet de ideale omgeving is om de modellering uit t e voeren. Een vertaalslag naar een programmeertaal als

C++

kan dit probleem verhelpen.

o Een tweede reden om over te gaan op de programmeertaal C++, is d a t de koppeling direkt binnen BARD (geschreven in

C++)

kan worden uitgevoerd. Dit zal leiden t o t een betere efficiëntie van de BEM en maakt het mogelijk om gebruik t e maken van slechts één invoer-file.

o Onderzocht moet worden of het gebaffelde cilinder model gebruikt kan worden voor de modellering van de akoestische afstraling van een MRI-systeem. Dit onderzoek vereist mogelijk de eerste prioriteit.

Bibliografie

[i] B. Soenarko. A boundary element formulation for radiation of acoustic waves from axi- symmetric bodies with arbitrary boundary conditions. Journal of the Acoustical Society of America, Vol. 93, Nr. 2, pag. 631-639,1993.

[2] A.D. Pierce. A C 0 USTICS A n Introduction to Its Physical Principles and Applications. McGraw-Hill Book Company, New York, 1981, ISBN 0-07-049961-6.

[3] C.A. Brebbia en J. Dominguez. Boundary Elements A n Introductory Course. Computa- tional Mechanics Publications, Southampton, 1992, ISBN 1-85312-160-6.

[4] R.D. Ciskowski en C.A. Brebbia. Boundary Element Methods in Acoustics. Elsevier Ap- plied Science, London, New York, 1991, ISBN 1085166-679-6.

[5] A.F. Seybert, C.Y.R. Cheng en T.W. Wu. The solution of coupled interior/exterior acoustic problems using the boundary element method. Journal of the Acoustical Society of America, Vol. 88, Nr. 3, pag. 1612-1618,1990.

Bijlage A

Opbouw

van

het BEM-pakket

(MATLAB-procedures)

De geschreven procedures kunnen grofweg in 7 groepen worden verdeeld: o hoofdprogramma (Rayleigh function)

o creëren van mesh-procedures (mesh function) e plot-procedures (plot function)

o vormfunctie-procedures (shape function) o hulp-procedures (utility function)

o koppel-procedures (utility function)

o extra-procedures (ten behoeve van controle) (extra function)

In figuur A.l staat het schema uit hoofdstuk 4 nogmaals weergegeven. Nu wordt echter aan- gegeven welke procedures achtereenvolgens worden aangeroepen wanneer het schema wordt doorlopen. Deze procedures zijn op diskette bijgevoegd. Binnen MATLAB kan met het commando: ‘HELP procedure-name’, een korte verklaring van de functie gevraagd worden. Als voorwaarden voor het kunnen werken met RAYLE1GH.M’ dient eerst het programma START.M te worden opgestart. Dit programma zorgt voor het actief maken van de juiste werk-directory’s. Vervolgens kan RAYLE1GH.M worden gestart. Bij wijziging van directory- structuur dienen zowel in START.M als in RAYLE1GH.M de juiste directory’s te worden gezet.

38 Opbouw van het BEM-pakket (MATLAB-procedures)

geheqgen.

ass.m

-

Bijlage

B

Numerieke integratie

B.1 Formule van Eder (trapeziumregel)

Interval ( a , b ) is verdeelt in n gelijke intervallen, zoals in figuur B.l staat weergegeven.

Figiiur B.l: Numeriek integreren

Bij een toename van het aantal intervallen neemt de nauwkeurigheid van de integratie toe.

B.2

Formule van Simpson

40 Numerieke integratie

Het aantal intervallen moet even zijn! Ook hier neemt de nauwkeurigheid van d e integratie toe bij een toenemend aantal intervallen.

B.3 QUAD8

-

standaard MATLAB-procedure

De procedure QUAD8.M is een hogere orde numerieke integratie methode. QUAD8 biedt de mogelijkheid om een tolerantie (standaard = 1

.

op t e geven. Deze procedure houdt zodoende rekening met het verloop van de functie.

B.4

GUASSHM

-

BARD procedure

Deze GAUSSHM.M routine berust op een Gauss-integratie (besproken in hoofdstuk 4). Ook hier geldt d a t bij het vergroten van het aantal intervallen de nauwkeurigheid toeneemt. De helling van de functie H , is daarbij van invloed op de nauwkeurigheid. Een toename van de helling leidt tot een afname van de nauwkeurigheid. Zodoende wordt binnen deze routine bij een stijler wordende helling het aantal intervallen vergroot. Het vergroten van het aantal intervallen is ook noodzakelijk, wanneer het golfgetal

IC

en de Fourier-coëfficiënt n toeneemt, omdat de periodetijd van H , dan kleiner wordt.

Bijlage C

Waarden van

en

wint

ten behoeve

van

Gauss-Legendre-integratie

n i p 1 - 2 3 - 4 5 - 6 - 7 .tip 0.00000 00000 00000 -0.57735 02691 89626 0.57735 02691 89626 -0.77459 66692 41483 0.00000 00000 00000 0.77459 66692 41483 -0.86113 63115 94053 -0.33998 10435 84856 0.33998 10435 84856 0.86113 63115 94053 -0.90617 98459 38664 -0.53846 93101 05683 0.00000 00000 00000 0.53846 93101 05683 0.90617 98459 38664 -0.93246 95142 03152 -0.23861 91860 83197 0.23861 91860 83197 0.66120 93864 66265 0.93246 95142 03152 -0.74153 11855 99394 -0.40584 51513 77397 0.00000 00000 00000 0.40584 51513 77397 0.74153 11855 99394 0.94910 79123 42759 -0.66120 93864 66265 -0.94910 79123 42759 - W i n t ,ip 2.00000 00000 00000 1.00000 00000 00000 1.00000 00000 00000 0.88888 88888 88888 0.55555 55555 55555 0.88888 88888 88888 0.34785 48451 37454 0.65214 51548 62546 0.65214 51548 62546 0.34785 48451 37454 0.23692 68850 56189 0.47862 86704 99366 0.55888 88888 88889 0.47862 86704 99366 0.23692 68850 56189 0.17132 44923 79170 0.36076 15730 48139 0.46791 39345 72691 0.46791 39345 72691 0.36076 15730 48139 0.17132 44923 79170 0.12948 49661 68870 0.27970 53914 89277 0.38183 00505 05119 0.41795 91836 73469 0.38183 00505 05119 0.27970 53914 89277 0.12948 49661 68870

42 Waarden van

E

en Utint ten behoeve van Gauss-Legendre-integratie .tip -0.96028 98564 97536 -0.79656 64774 13627 -0.52553 24099 16329 -0.18343 46424 95650 0.18343 46424 95650 0.52553 24099 16329 0.79666 64774 13627 0.96028 98564 97536 -0.96816 02395 07626 -0.83603 11073 26636 -0.61337 14327 00590 -0.32425 34234 03809 0.00000 00000 00000 0.32425 34234 03809 0.61337 14327 00590 0.83603 11073 26636 0.96816 02395 07626 -0.97390 65285 17172 -0.86506 33666 88985 -0.67940 95682 99024 -0.43339 53941 29247 -0.14887 43389 81631 0.14887 43389 81631 0.43339 53941 29247 0.67940 95682 99024 0.86506 33666 88985 0.97390 65285 17172 -0.98156 06342 46719 -0.90411 72563 70475 -0.?6990 26741 94305 -0.58731 79542 86617 -0.36783 14989 98180 0.12523 34085 11469 0.36783 14989 98180 0.58731 79542 86617 0.76990 26741 94305 0.90411 72563 70475 0.98156 06342 46719 -0.12523 34085 11469 .Wi,t,ip 0.10122 85362 90376 0.22238 10344 53374 o.313'70 66458 7ï887 0.36268 37833 78362 0.36268 37833 78362 0.31370 66458 77887 0.22238 10344 53374 0.10122 85362 90376 0.08127 43883 61574 0.18064 81606 94857 0.26061 06964 02935 0.31234 70770 40003 0.33023 93550 01260 0.31234 70770 40003 0.26061 06964 02935 0.18064 81606 94857 0.08127 43883 61574 0.06667 13443 08688 0.14945 13491 50581 0.21908 63625 15982 0.26926 67193 09996 0.29552 42247 14753 0.29552 42247 14753 0.26926 67193 09996 0.21908 63625 15982 0.14945 13491 50581 0.06667 13443 08688 0.04717 53363 86512 0.10693 93259 95318 0.16007 83285 43346 0.20316 74267 23066 0.23349 25365 38355 0.24914 70458 13403 0.24914 70458 13403 0.23349 25365 38355 0.20316 74267 23066 0.16007 83285 43346 0.10693 93259 95318 0.04717 53363 86512 Tabel C.l: Gauss-Legendre-integratie

Bijlage

D

Log-Gauss-transformatie

Log-Gauss-integratie is van de vorm:

Deze integratie methode willen we transformeren tot in de vorm van de Rayleigh integraal.

D . l

Geval

1

In figuur D.l staan de beide locale coördinatensystemen van geval 1 weergegeven. Met behulp

geval 1

-

5 - 1 5=0 g=i

q = O q=l

@=knooppunt )<=singulariteit

Figuur D.l: Geval 1

hiervan kunnen de volgende transformatie-vergelijkingen worden afgeleid. Tevens is ook diens afgeleide bepaald:

D.3 Geval 3 45

D.3

Geval 3

In figuur D.3 staan de beide locale coördinatensystemen van geval 3 weergegeven. Nu gelden

geval 3

i

5 - 1

5=0

g=i q=1 17=0 .=knooppunt X=singulariteit Figuur D.3: Geval 3 de volgende transformatie-vergelijkingen en afgeleide:

q = - ( - [ + l )

1

¢j t = - 2 7 + 1 ,

2

Substitutie van deze uitdrukkingen in de definitie vergelijking van Log-Gauss levert ( f ( q i p ) =

s ( - 2 a e p + l )

.

- w ï Ï x ) *

Waarden van

q en

wint

ten behoeve

van Log-Gauss-integratie

0.60227691 0.36899706 0.76688030 ~ 0.41448480e-1

:0:;;::;545

I 0.84898239 O 2 9 134472e-1 0.17397721

,

0.41170251 O . 6 ? 7 3 M 7 0.89477136 0.216344005e-1 0.12958339 0.31402045 0.53865721 0.75691533 0.92266884 0.16719355e-1 0.10018568 0.24629424 0.43346349 0.63235098 0.81 11 1862 O .940848 16 -Wint,ip 0.71853931 0.28146068 0.51340455 0.39198004 0.94615406e-1 0.38346406 0.38687532 0.19043513 0.39225487e-1 0.29789346 0.34977622 0.23448829 O .98!XWO6Oe- 1 0.18911552e-1 0.23876366 0.30828657 0.24531742 O. 14200875 O .55454622e- 1 0.10168958e-1 0.19616938 0.27030264 0.23968187 O. 16577577 0.88943226e-1 O .33 194304e- 1 0.59327869e-2