Zweidimensionales, akustisches Design zwischen Schalenelementen eines

Flugzeugrumpfs mittels Topologieoptimierung

Martin Radestock1, Michael Rose2, Hans Peter Monner31 Deutsches Zentrum für Luft- und Raumfahrt, 38108 Braunschweig, Deutschland, Email: Martin.Radestock@dlr.de

2 Deutsches Zentrum für Luft- und Raumfahrt, 38108 Braunschweig, Deutschland, Email: Michael.Rose@dlr.de3 Deutsches Zentrum für Luft- und Raumfahrt, 38108 Braunschweig, Deutschland, Email: Hans.Monner@dlr.de

EinleitungBei der Findung von optimalen Leichtbaukonstruktio-nen kommt es oftmals auf die Erfahrung von beteilig-ten Personen an. Lassen sich diese Problemstellungenmit Ihren entsprechenden Randbedingungen mathema-tisch vollständig beschreiben, so kann mit Hilfe eines Al-gorithmus ein optimales Design für ein Bauteil gefun-den werden, das bei minimalem Materialeinsatz den äu-ßeren Belastungen standhält. Bei diesem Verfahren han-delt es sich um die Topologieoptimierung, welche sehrhäufig im Rahmen von Strukturoptimierungen eingesetztwird[1],[2]. Da es einen direkten Zusammenhang zwi-schen der Elementsteifigkeit und dem verteilten Materialgibt, ist dieses Verfahren im Rahmen der Steifigkeitsop-timierung sehr ausgereift und weit verbreitet. Auf Grunddes bekanntermaßen häufig ungünstigen akustischen Ver-haltens von Leichtbaukonstruktion, ist eine Anpassungder Zielfunktion auf akustische Problemstellungen vonhohem Interesse. Im Rahmen dieses Artikels wird die-se Anpassung beschrieben. Dazu zählt die Beschreibungdes Materials im Zusammenhang mit der Helmholtz-Differentialgleichung und die Formulierung der Topolo-gieoptimierung. Abschließend werden zwei Beispiele vor-gestellt und die Ergebnisse diskutiert.

Grundlegende akustische Annahmen fürdie OptimierungDer Ansatz für die Akustik im zweidimensionalen Raumist die Helmholtz-Differentialgleichung [3].

∆p(r) + k2p(r) = 0 (1)

Hierbei ist ∆ der Laplace Operator, p der Schalldruck, kdie Wellenzahl und r die Ortskoordinate für den Schall-druck. Für die Optimierung wird der Schalldruck als zuminimierende Schallgröße gewählt. Die Wellenzahl k ent-hält die Kennkreisfrequenz ω und die Schallgeschwindig-keit c, welche vom Material abhängt.

k =ω

c

[1

m

](2)

Für die Materialverteilung muss die Schallgeschwindig-keit für Luft und ein festes Material beschrieben werden.Die Schallgeschwindigkeit in Luft cgas ist über den adia-batischen Kompressionsmodul K und die Dichte ρ defi-niert.

cgas =

√K

ρ

[ms

](3)

Bei dem festen Material muss die Schallgeschwindigkeitauch nur durch Materialgrößen beschreibbar sein. Dasgilt jedoch nur für die Longitudinalgeschwindigkeit csolid,die über das Elastizitätsmodul E, die Dichte ρ und dieQuerkontraktionszahl ν definiert ist.

csolid =

√E(1− ν)

ρ(1 + ν)(1− 2ν)

[ms

](4)

Die verwendeten Materialien sind zum einen Luft mitder Dichte ρLuft = 1.2 kg/m3, dem KompressionsmodulKLuft = 1.42e5 Pa und als Festkörper Aluminium mitder Dichte ρAl = 2700 kg/m3 , dem ElastizitätsmodulE = 70e9 Pa und der Querkontraktionszahl ν = 0.33. Ei-ne direkte Beziehung zwischen Kompressionsmodul undElastizitätsmodul bzw. zwischen den Dichten in den Glei-chungen (3) und (4) kann nicht hergestellt werden. Durchdie Vernachlässigung der Querkontraktionszahl entstehtein Fehler von ca. 16 % bei der Schallgeschwindigkeit imAluminium. Die Designvariable ξ beschreibt in der Opti-mierung den Übergang zwischen der Schallgeschwindig-keit von Luft und Aluminium. In einem Wertebereich von0 bis 1 wird für die Designvariable folgende Beziehungverwendet.

c = cgas + (csolid − cgas)ξ(r)p

[ms

](5)

Der Faktor p ist ein Bestrafungsfaktor, der die Schallge-schwindigkeit stärker zu Ihren Randwerten zwingt. DieDesignvariable variiert kontinuierlich in ihrem Wertebe-reich und durch den Bestrafungsfaktor entstehen weni-ger Zonen mit Zwischenmaterial, deren Schallgeschwin-digkeit zwischen der von Luft und Aluminium liegt [4].

Formulierung der TopologieoptimierungIn mathematische Form wird die Topologieoptimierungmit Zielfunktion und Nebenbedingungen wie folgt zusam-mengefasst.

minξ δ = minξ1∫

Ξtardr

∫Ξtar

|p(r, ξ(r))|2dr

(6)

s.t.1∫

Ξdesdr

∫Ξdes

ξ(r)

βdr − 1 ≤ 0 (7)

0 ≤ ξ(r) ≤ 1 (8)

DAGA 2017 Kiel

729

Hierbei ist δ die Zielfunktion, Ξtar der Zielbereich aufden ein Raum optimiert wird, Ξdes ist der Designbereichin dem Material verteilt werden kann und β der Füll-grad für den Designbereich. Die Formulierung ist ähn-lich zur bestehenden Literatur [5], [6]. Es wird als Ziel-funktion nur lokal der Schalldruck minimiert und nichtglobal in der gesamten Raum. Am Beispiel eines einfa-chen Rechteckraums werden die relevanten Bereiche ausder mathematischen Beschreibung erklärt. Der Raum istin Abbildung 1 dargestellt. Der gelb umrandete Bereichstellt mit den Elementen den zu minimierenden Zielbe-reich Ξtar dar. Im oberen Bereich dieses Raumes ist eindunkelgrauer Bereich Ξdes definiert, in dem die Mate-rialverteilung erfolgen soll. Nur diese Elemente könnenkontinuierlich verändert werden. Alle anderen Elementebesitzen die Eigenschaft von Luft. Die Erregung für die-sen Raum wird durch eine Punktquelle mit dem Radiusvon 0,1m und einer Oberflächengeschwindigkeit von 0,01m/s modelliert.

Abbildung 1: Vernetzter 2D Raum für die Topologieopti-mierung

Das Netz wird mit Hilfe einer Delauny Triangulation be-rechnet [7]. Hierbei können u.a. bestimmte Randbereiche,der maximale Flächeninhalt und der minimale Dreiecks-winkel vorgegeben werden. Weiterhin besitzen die Ele-mente eine quadratische Ansatzfunktion mit der die FEBerechnung stattfindet. Zur Aktualisieren der Elemen-te wird ein gängiger Algorithmus in der Topologieopti-mierung verwendet. Die ”method of moving asymptotes”(MMA) ist ein gradientenbasiertes Verfahren, welches dieZielfunktion und die Nebenbedingungen berücksichtigt[8].

Sensitivitätsanalyse

Für die Aktualisierung der Variablen ist aufgrund desgradientenbasierten Verfahrens eine Sensitivität notwen-dig. Diese muss sowohl für die Zielfunktion als auch fürdie Nebenbedingungen erstellt werden. Allgemein ist dieBeschreibung für die FE Berechnung linear und kannüber eine komplexe Systemmatrix S, den komplexenDruckvektor p und dem komplexen Anregungsvektor fausgedrückt werden. Die fett gekennzeichneten Symbole

markieren Vektoren oder Matrizen.

S(ξ)p(ξ) = f (9)

Aufgrund der impliziten Beziehung zwischen Designva-riable und Zielfunktion wird die adjoint state methodangewendet [9]. Dazu wird die Zielfunktion mit Nullter-men erweitert.

δ = δ + λT (Sp− f) + λT(Sp− f) (10)

Hierbei ist λ ein komplexer Vektor, der zunächst freiwählbar ist, T markiert den transponierten Vektor unddie überstrichene Variablen sind komplex konjugierteGrößen. Die Sensitivität für die Zielfunktion δ lässt sichin der folgenden Form beschreiben. Die Berechnung er-folgt in Anlehnung an die Literatur [5],[6].

dδ

dξ=

∂δ

∂ξ+

∂δ

∂pre

∂pre

∂ξ+

∂δ

∂pim

∂pim

∂ξ(11)

Der Anteil ”re” beschreibt den Realteil des komplexenDruckvektors und ”im” den Imaginärteil. Um partiel-len Ableitungen zu entgehen, die nicht bestimmbar sind,wird der komplexe Vektor λ wie folgt gewählt.

λTS = −

(1∫

Ξtardr

∫Ξtar

(pre − jpim)dr

)(12)

Um den komplexen Vektor zu ermitteln, wird eine weitereRechnung durchgeführt. Vereinfacht ausgedrückt wird imZielbereich Ξtar die negativ komplex konjugierte Lösungder Ausgangsrechnung als Anregungsvektor genutzt, umeine zweite Lösung für das Schalldruckfeld zu generieren.Die Sensitivität vereinfacht sich mit dem Vektor λ dannzur folgenden Formulierung.

dδ

dξ=

∂δ

∂ξ+ 2Re

[λT ∂S

∂ξp

](13)

Auch zur Nebenbedingung muss eine Ableitung gebildetwerden. Diese lässt sich direkt beschreiben und hängt fürjedes Element im Designbereich Ξdes von der Volumen-grenze ab.

d

dξ

(1∫

Ξdesdr

∫Ξdes

ξ(r)

βdr − 1

)=

1

β(14)

FE AnalyseEine wesentlich Herausforderung für die Sensitivitätsana-lyse ist die partielle Ableitung der Systemmatrix nachder Designvariablen. Da die meisten kommerziellen Pro-gramme, wie z.B. ANSYS®, über diese Funktion nichtverfügen, muss eine eigene FE-Berechnung implemen-tiert werden. Dazu wurde in einer MATLAB®-Umgebungdie Beschreibung objektorientiert über ein MEX-plugindirekt bereit gestellt. Dieses Vorgehen wurde gewählt,um die MATLAB® Umgebung mit ihren Vorteilen inVektor- und Matrixoperationen auszunutzen. Die ma-thematische Beschreibung zur Umsetzung der Helholtz-Differentialgleichung für ein FE-Netz folgt der Literatur[10] und wird nicht weiter beschrieben.

DAGA 2017 Kiel

730

SimulationsbeispieleDas erste Simulationsbeispiel ist der Rechteckraum ausAbbildung 1. Dazu sind zum einen der Amplituden-Frequenzgang in Abbildung 2 und die Eigenfrequenz desRaumes bei 34.45 Hz in Abbildung 3 dargestellt. DieseFrequenz wurde frei gewählt und es soll gezeigt werden,dass der Schalldruck in dem Raum bei dieser Eigenfre-quenz für den eingekreisten Bereich reduziert wird.

Abbildung 2: Amplituden-Frequenzgang des Rechteckrau-mes ohne eine Materialverteilung

Abbildung 3: Eigenfrequenz des Rechteckraumes bei 34.45Hz ohne Materialverteilung

Der Amplituden-Frequenzgang für das optimierte De-sign ist in Abbildung 4 und die Schalldruckverteilungbei 34.45 Hz ist in Abbildung 5 dargestellt. Es zeigensich klar definierte, schwarze Elemente, in denen Ma-terial vorhanden ist. Hier zeigt sich jedoch schon eineHerausforderung, wenn diese Strukturen gefertigt wer-den sollen. Es gibt Elemente, die frei im Raum liegenund keine Verbindung zu anderen Elementen besitzen.Außerdem lässt sich im Amplituden-Frequenzspektrumerkennen, dass die Eigenfrequenz sich durch die Materi-alverteilung verschoben hat.

Abbildung 4: Amplituden-Frequenzgang des Rechteckrau-mes mit optimierter Materialverteilung

Daher ist die Anwendung der Topologieoptimierung eherfür tonale Anwendungen oder für ein schmales Frequenz-band zu sehen, in dem Spitzen im Schalldruck reduziertwerden können. Das zweite Beispiel bezieht sich daher

Abbildung 5: Schalldruckverteilung des Rechteckraums bei34.45 Hz mit optimierter Materialverteilung

auf einen Flugzeugrumpf ähnlich einer Airbus 320, dietonal angeregt werden kann. (Abbildung 6). Die Erre-gung soll der power control unit(PCU) entsprechen, diezur Steuerung der Klappen und Steuerflächen eine tonaleAnregung bei ca. 800 Hz erzeugt. Hier soll im approxi-

Abbildung 6: Aufbau der Rumpftonne

mierten Kopfbereich der Passagiere der Schalldruck ver-mindert werden. Dazu kann im rot abgegrenzten BereichMaterial verteilt werden, in dem das Kabinenlining sichbefindet. Die blauen Begrenzungen zeigen die Rumpfton-ne und den Boden des Flugzeugs und werden im Rahmender Topologieoptimierung konstant als volles Material ge-rechnet. Die Erregung wird durch eine konstante Ober-flächengeschwindigkeit realisiert. Diese ist nicht identischmit der Anregungsstärke durch die PCU. In Abbildung7 ist der mittlere Schalldruck im Zielbereich dargestellt,in dem sichdrei Spitzen ausbilden. Der Schalldruck bei800 Hz ist in Abbildung 8 ohne Materialverteilung dar-gestellt.

Abbildung 7: Amplituden-Frequenzgang des mittlerenSchalldrucks im Zielbereich der Rumpftonne

Für die Optimierung wurde der Frequenzbereich von

DAGA 2017 Kiel

731

Abbildung 8: Schalldruckverteilung in der Tonne bei 800 Hz

797 Hz bis 803 Hz vorgegeben, da hier die Spitzenim Amplituden-Frequenzgang minimiert werden sollen.Abbildung 9 zeigt sowohl die vorgegebenen Frequenzenals auch den Amplituden-Frequenzgang der optimiertenStruktur. Wichtig ist daher auch die Schalldruckvertei-lung für die maximal auftretende Amplitude, die in Ab-bildung 10 dargestellt ist.

Abbildung 9: Amplituden-Frequenzgang des mittlerenSchalldrucks im Zielbereich mit optimierter Materialvertei-lung

Abbildung 10: Schalldruckverteilung in der Tonne für diemaximal auftretende Amplitude bei 802 Hz

Das Ergebnis ist bei diesem Beispiel nicht so eindeutigwie zuvor. Es gibt noch graue Bereiche und auch losgelö-ste Teile, die frei im Raum hängen. Jedoch zeigt die Ma-terialverteilung schon den Einfluss auf den Schalldruck,wobei hier nur die Helholtz-Differentialgleichung umge-setzt wurde. Effekte durch unterschiedliche Impedanzenwurden bisher noch nicht berücksichtigt. Der Übergangan Elementen mit unterschiedlichen Impedanzen verur-

sacht u.a. eine teilweise Reflektion, die durch reine Va-riation der Schallgeschwindigkeiten in der Helmholtz-Differentialgleichung nicht abgedeckt ist.

FazitAuf Basis einer Topologieoptimierung wurde eine akusti-sche Optimierung erzielt, die mit einer MaterialverteilungEinfluss auf den Schalldruck eines Raumes nimmt. Da-bei wird direkt die Helmholtz-Differentialgleichung ver-ändert, wobei die Materialverteilung voranging die Ver-schiebung von z.B. Eigenfrequenzen eines geschlossenenRaumes vornimmt. Am Beispiel eines Rechteckraumsund eines approximierten Flugzeugrumpfes mit tona-ler Erregung wird die Topologieoptimierung angewendet.Weiterführende Arbeiten sind die Umsetzung von Im-pedanzübergängen in der Optimierung und Algorithmenzur Vermeidung von frei liegendem Material.

Literatur[1] Zhu, J.-H., Zhang, W.-H. und Xia, L.: Topology opti-

mization in aircraft and aerospace structures design.Archives of Computational Methods in Engineering(2016), 595–622

[2] Krog, L., Tucker, A., Kemp, M. und Boyd, R.: To-pology optimisation of aircraft wing box ribs. 10thAIAA/ISSMO multidisciplinary analysis and opti-mization conference (2004), 4481 .

[3] Kuttruff, H.: Room acoustics. Spon Press, London,2009

[4] Bendsøe, M. P. und Sigmund, O.: Topology optimiza-tion: Theory, methods, and applications. Springer,Berlin, 2004

[5] Dühring, M. B., Jensen, J. S. und Sigmund, O.:Acoustic design by topology optimization. Journal ofSound and Vibration (2008), 557–575

[6] Jensen, J. S. and Sigmund, O.: Systematic design ofacoustic devices by topology optimization. ICSV12Proceedings (2005)

[7] Shewchuk, J. R.: Triangle: Engineering a 2D QualityMesh Generator and Delaunay Triangulator. AppliedComputational Geometry: Towards Geometric Engi-neering (1996), 203–222, URL:http://www.cs.cmu.edu/~quake/triangle.html

[8] Svanberg, K.: The method of moving asymptotes—anew method for structural optimization. Internatio-nal journal for numerical methods in engineering(1987), 359–373

[9] Tortorelli, D. A. und Michaleris, P.: Design sensitivityanalysis: Overview and review. Inverse Problems inEngineering (1994), 71–105

[10] Ihlenburg, F.(Ed.): Finite Element Analysis of Acou-stic Scattering. Springer-Verlag New York Inc, NewYork, 1998

DAGA 2017 Kiel

732