FDS: Numerische Optimierung und Kooperationen · Fassung 21. September 2009 autoren DR. RER. NAT....

Ingenieure für Brandschutz

FDS: Numerische Optimierung

und Kooperationen

Dr. rer. nat. Susanne Kilian

Dr.-Ing. Christoph Klinzmann

hhpberlin – Ingenieure für Brandschutz GmbH

Fassung September 2011

InHalt

1 Einleitung

2 Zielsetzungen und Aktivitäten der FDS Usergroup

3 Optimierte Parallelisierungsstrategien für den FDS-Drucklöser

3.1 Die FDS-Druckgleichung

3.2 Algorithmische Beschreibung des aktuellen FFT-Drucklösers

3.3 Algorithmische Beschreibung des alternativen

ScaRC-Drucklösers

3.4 Entwicklungsstand des ScaRC-Drucklösers

4 Numerische Testrechnungen

4.1 Durchströmung eines Kanals in 2D

4.2 Dreifach-Einströmung in ein kubisches Gebiet in 3D

4.3 Diskussion der Testreihen

5 Fazit

Fassung 21. September 2009

autoren

DR. RER. NAT. SUSANNE KIlIAN

Dr. Susanne Kilian ist Mathematikerin mit Schwerpunkt Numerik und

kann auf eine langjährige akademische Erfahrung in der Entwicklung

paralleler Programmpakete zur Simulation strömungsmechanischer

Probleme (CFD) auf Hochleistungscomputern zurückgreifen. Bei

hhpberlin fördert Sie die Entwicklung und Anwendung von FDS für

den schutzzielorientierten Brandschutz. Gemeinsam mit hhpberlin ini-

tiierte Dr. Kilian zu diesem Zweck 2008 die FDS Usergroup und forscht

aktiv an der Weiterentlicklung von FDS.

DR.-ING. CHRISTOPH KlINZMANN

Dr.-Ing. Christoph Klinzmann unterstützt von unserem Hamburger

Standort aus bundesweite Bauvorhaben mittels wissenschaft-

lichen Untersuchungen und Anwendungen. Zudem betreut er die

Konzipierung und Umsetzung von Brandschutzmaßnahmen u. a.

für Bauprojekte von Airbus und ECE. Im Normenausschuss DIN

18230-4 „Brandsimulation“ trägt er dazu bei, den Brandschutz

weiterzuentwickeln und zu optimieren.

6


FDS: nuMerISCHe oPtIMIerunG

unD KooPeratIonen

1 EINlEITUNG

Immer komplexer werdende Gebäude mit architektonischen Besonderheiten,

wie etwa Atrien, stellen Ingenieure und Bauaufsichtsbehörden häufig vor die

Frage, wie die zugrunde liegenden Gebäudestrukturen und ihre Abweichungen

vom üblichen Standardgebäude mit dem starren Baurecht in Einklang gebracht

werden können. Speziell im Fall der ausreichenden Rauchableitung oder der

Standsicherheit im Brandfall muss bei diesen Gebäuden der Nachweis über

moderne Ingenieurmethoden erfolgen.

Computergestützte Simulationsprogramme auf Basis moderner Techniken der

Computational Fluid Dynamics (CFD) ermöglichen heutzutage in zunehmendem

Maße eine individuelle und maßgeschneiderte Abbildung der realen Gebäude-

abmessungen und der maßgeblichen Brandszenarien. Ein Beispiel für ein re-

nommiertes CFD-Programm zur numerischen Brand- und Rauchgas-Simulation,

das weltweit auf breiter Basis eingesetzt wird, ist der Fire Dynamics Simulator

(FDS) [1]. FDS wurde am National Institute of Standards and Technology (NIST)

in Maryland/USA entwickelt und steht unter einer Open-Source lizenz zur

freien Verfügung.

Üblicherweise handelt es sich bei den genannten CFD-Nachweisen um aus-

gesprochen laufzeit- und Speicherplatz-intensive Simulationen, die nur auf

modernen Parallelrechner-Architekturen effizient durchgeführt werden können.

Die aktuell eingesetzten CFD-Programme basieren dabei auf hochkomplexen,

numerischen Methoden, die über viele Jahre oder gar Jahrzehnte hinweg entwi-

ckelt wurden. Der eingesetzte Programmcode entstammt zumeist aus Zeiten, in

denen Parallelrechner-Architekturen noch nicht weit verbreitet waren. So liegen

die ausgeklügelten und über einen langen Zeitraum bewährten numerischen

Methoden zur Berechnung der betrachteten physikalischen Größen vielfach nur

für Einzelprozessor-Anwendungen vor. Die numerische Anpassung an die neuen

Multiprozessor-Architekturen ist ein höchst sensibler Prozess, der einen großen

strukturellen Aufwand erfordert. Insbesondere muss streng darauf geachtet

werden, dass die Konvergenzeigenschaften und die Approximationsgenauigkeit

der ursprünglichen seriellen Methodik nicht beeinträchtigt werden.

7

FDS: Numerische Optimierung und Kooperationen

In diesem Zusammenhang stellt insbesondere die effiziente Parallelisierung der

FDS-Druckgleichung noch eine große Herausforderung dar. Um Genauigkeits-

verluste bzw. numerische Instabilitäten zu vermeiden, muss die verwendete

Strategie zur Parallelisierung unbedingt den ausgeprägt globalen Charakter

des physikalischen Drucks berücksichtigen.

Die gegenwärtig in FDS eingesetzte Methodik zur lösung der Druckgleichung,

die parallele Durchführung verteilter Fast Fourier Transformationen (FFT),

weist einen eher grobkörnigen und lokalen Charakter auf und wird dieser

Anforderung leider nicht in allen Fällen gerecht. So kann insbesondere die

fehlerfreie Skalierung auf große Teilgebietszahlen nicht zuverlässig gewähr-

leistet werden. Die Verbesserung der Qualität der parallelen FDS-Version,

speziell im Hinblick auf massiv parallele Anwendungen mit hunderten von

Prozessoren, wird daher in der offiziellen FDS-Roadmap [2] als sehr dring-

liches Ziel beschrieben.

In enger Kooperation mit dem amerikanischen FDS-Entwicklerteam am NIST

widmet sich hhpberlin daher seit mehreren Jahren dem Design und der

Implementierung des alternativen parallelen Drucklösers Scalable Recursive

Clustering (ScaRC). Dieser basiert im Gegensatz zum aktuellen FFT-Drucklöser

auf einer globalen Philosophie, wodurch gebietsübergreifende Problemabhän-

gigkeiten bestmöglich reproduziert werden.

Begleitet wird diese enge Zusammenarbeit durch die deutsche FDS Usergroup,

die als übergreifende Plattform für einen intensiven anwenderorientierten und

wissenschaftlichen Austausch rund um FDS ins leben gerufen wurde. Eine

wichtige Zielsetzung der FDS Usergroup bestand von Anfang an darin, einen

persönlichen Kontakt zwischen den Mitgliedern der FDS Usergroup und den

FDS-Entwicklern herzustellen. Auf diese Weise war es in jüngster Vergangen-

heit bereits mehrfach möglich, Fragen und Anregungen direkt an die Haupt-

entwickler weiterzugeben und so Einfluss auf die laufende Entwicklung von

FDS zu nehmen.

Dank des kompetenten und engagierten Entwicklerteams rund um Kevin

McGrattan (NIST) ist FDS einer stetigen Fortentwicklung unterworfen. In der

offiziellen FDS-Roadmap [2] haben die Hauptentwickler eine Reihe kurz- und

mittelfristiger Ziele definiert, zu deren Umsetzung die Unterstützung durch

interessierte Entwickler-Teams aus aller Welt sehr willkommen ist.

8


Der vorliegende Artikel gliedert sich wie folgt: Das erste Kapitel gibt zunächst

eine kurze Zusammenfassung über die Ziele und Aktivitäten der deutschen

FDS Usergroup. Um einen Überblick über die grundlegenden Funktionsweisen

des FFT- bzw. ScaRC-Zugangs zu vermitteln, schließt sich daran eine informelle

algorithmische Beschreibung beider Verfahrensklassen an. Es folgen eine

Darstellung des aktuellen ScaRC-Entwicklungsstands sowie ein Ausblick auf

die in naher Zukunft geplanten Entwicklungsschritte. Um den Zugewinn an

Approximationsgenauigkeit und Skalierbarkeit durch die neue ScaRC-Methodik

zu illustrieren, werden nachfolgend zwei numerische Vergleichsrechnungen prä-

sentiert, die die betrachteten Drucklöser einander gegenüberstellen. Der Artikel

schließt mit einer Diskussion der numerischen Resultate und einer Bewertung

der bisher erzielten Erfolge.

2 ZIElSETZUNGEN UND AKTIVITÄTEN DER FDS USERGROUP

Die FDS Usergroup wurde im November 2007 auf die Initiative von hhpberlin

hin gegründet und ist seitdem stetig angewachsen [4]. Mittlerweile sind zahl-

reiche Ingenieurbüros und Forschungseinrichtungen Mitglied dieser Gruppe.

Die FDS Usergroup orientiert sich an den folgenden Zielsetzungen:

• aufbau eines internationalen Kontakt- und Kooperationsnetzwerks:

Als unabhängiges Forum für FDS Anwender und Forscher aus deutschspra-

chigen und zunehmend auch anderen europäischen ländern versteht sich

die FDS Usergroup als internationale Plattform zum wissenschaftlichen

Austausch von Informationen, Erfahrungen und Fragen rund um FDS.

Unter dem Motto Share Know-how sollen Talente aus allen Bereichen der

Anwendung, Verifikation, Validierung und numerischen Weiterentwicklung

von FDS miteinander verknüpft werden, um gemeinsam zu einer Qualitäts-

steigerung von FDS beizutragen.

• Koordination der FDS-aktivitäten: Eine besondere Bedeutung hat die

gemeinsame Identifikation des weiteren Forschungs- und Entwicklungs-

bedarfs sowie die Definition einheitlicher Qualitätsstandards und Anwen-

dungsrichtlinien im Rahmen eines fortlaufenden Diskussionsprozesses.

• organisation von Workshops: Die Veranstaltung regelmäßiger Treffen und

Workshops soll den Aufbau und die Pflege persönlicher Kontakte zwischen

9


den einzelnen Mitgliedern erleichtern. Gleichzeitig sollen neueste Erkennt-

nisse einem breitem Publikum zugänglich und für die tägliche Anwen-

dungspraxis nutzbar gemacht werden.

• regelmäßiger Kontakt zu den FDS-Hauptentwicklern: Eine wichtige Ver-

mittlerfunktion besteht darin, die Interessen und Konzepte der Mitglieder

in gebündelter Weise an die FDS-Hauptentwickler heranzutragen, um ihnen

auf diese Weise mehr Nachdruck zu verleihen. Dies soll letztlich ermögli-

chen, bei der weiteren Entwicklung von FDS aktiv mitgestalten zu können.

Ursprünglich nur gedacht für den deutschsprachigen Raum, verzeichnet die

FDS Usergroup gegenwärtig 260 Mitglieder aus Deutschland, Österreich, der

Schweiz, Belgien, Frankreich, Norwegen, Dänemark und Finnland. Das Online-Fo-

rum, das registrierten Mitgliedern zur Verfügung steht, ermöglicht es, konkrete

Anwenderprobleme zeitnah zu diskutieren und gemeinsam zu lösen bzw. sich

über die Aktivitäten anderer Teilnehmer zu informieren.

Im Verlauf der zurück liegenden drei Jahre hat sich ein stabiler Kreis an Veran-

staltungs-Teilnehmern heraus gebildet, in dem man sich wechselseitig kennt und

voneinander profitiert. Wir blicken inzwischen zurück auf vier sehr erfolgreiche

Anwendertreffen sowie regelmäßige regionale Treffen und einen Workshop der

Arbeitsgruppe NRW, in deren Verlauf unterschiedlichste Themen rund um FDS

präsentiert und diskutiert wurden. Dabei reichte die Palette der dargestellten

Themenbereiche von Beispielen und Problemen aus der täglichen Ingenieurpra-

xis über die Darstellung neuer Parallelisierungs- und Verifikationskonzepte bis

hin zur Präsentation vollständig neuer Code-Bestandteile. Das fünfte Anwen-

dertreffen ist für November 2011 bereits in Planung. Dabei sind neue Mitglieder

jederzeit herzlich willkommen.

Die während des jüngsten Anwendertreffens im November 2010 eingerichteten

Arbeitskreise boten erstmalig die Möglichkeit, in themenorientierten Klein-

gruppen erste Konzepte für die gemeinsame Definition qualitätssichernder

Kriterien zu erarbeiten. Eine wichtige Zielsetzung, die sich aus dem Arbeitskreis

zur Qualitätssicherung ergeben hat, ist die Konzeption eines gemeinsamen

leitfadens, der die Grundlage für ein einheitliches Prozedere bei der Erstellung

und Dokumentation von FDS-Gutachten bilden soll. Ein grundsätzlich standar-

disiertes Vorgehen gegenüber Bauaufsichtsbehörden auf Basis eines gemein-

samen leitfadens kann ein wichtiges Instrument darstellen, um die Qualität

sowie Glaubwürdigkeit von FDS-Simulationen zu erhöhen. Eine erste Version

10


des leitfadens soll im Verlauf dieses Jahres konkrete Gestalt annehmen und auf

den kommenden Treffen ausführlich diskutiert werden.

Insbesondere die Roundtable Discussion im Juni 2010 hat neue Impulse gesetzt

und einen „heißen Draht“ zwischen einzelnen Mitgliedern der FDS Usergroup

ins leben gerufen: Die beiden FDS-Entwickler Kevin McGrattan und Randall

McDermott waren zu Gast in Berlin und boten allen Teilnehmern im Rahmen

einer offenen Diskussionsrunde die Möglichkeit, konkrete Fragen und Probleme

intensiv zu diskutieren und zu klären. Wie im Nachhinein von den Hauptentwick-

lern betont wurde, hat dieses Treffen unmittelbaren Einfluss auf die aktuelle

Weiterentwicklung von FDS genommen. Dies geht deutlich aus dem Grußwort

der Hauptentwickler Kevin McGrattan und Randall McDermott anlässlich des

vierten Anwendertreffens der FDS Usergroup hervor:

„ ... As you can see, several of these enhancements are a direct result of our

June meeting in Berlin and we are thankful for our wonderful collaboration ...“.

Diese intensive Zusammenarbeit soll auch zukünftig fortgeführt und durch

weitere regelmäßige Treffen untermauert werden.

3 OPTIMIERTE PARAllElISIERUNGSSTRATEGIEN FÜR DEN FDS-DRUCKlÖSER

Vor dem Hintergrund der rasanten Entwicklungen im Bereich der Multiprozessor-

Architekturen bietet FDS die Möglichkeit, Berechnungen parallel auf mehreren

Prozessoren auszuführen, wodurch potentiell eine deutliche Reduktion der

Gesamtrechenzeit erreicht werden kann. Zu diesem Zweck wird das zugrunde

liegende Berechnungsgebiet vom Anwender in kubische Teilgebiete mit recht-

eckigen Gittern unterteilt, die dann weitgehend unabhängig von den Prozes-

soren eines Parallelrechners bearbeitet werden können.

Da aufgrund von komplexen Strömungs- und Austauschprozessen die einzelnen

Teilgebiete nicht unabhängig voneinander betrachtet werden können, müssen

die zugehörigen Teilprozesse regelmäßig miteinander gekoppelt werden. Dies

erfolgt durch den koordinierten Austausch wechselseitig benötigter Daten auf

Basis der Kommunikationsbibliothek Message Passing Interface (MPI) [6]. Diese

Vorgehensweise setzt notwendigerweise einen stabilen parallelen lösungs-

11


algorithmus voraus, der dazu in der lage ist, die betrachtete Gebietszerlegung

zuverlässig und numerisch stabil zu bearbeiten.

Ein genauer Blick in die numerischen Grundlagen von FDS offenbart, dass der

zeitliche Verlauf der betrachteten thermodynamischen Größen unter Ver-

wendung expliziter Zeitschrittmethoden approximiert wird, siehe [15, 16]. Dies

bedeutet konkret, dass bei der Berechnung neuer Näherungswerte nur auf

zuvor bereits berechnete Werte aus vorangehenden Zeitschritten zurückge-

griffen wird. Da explizite Methoden nur einen lokalen Datenaustausch zwischen

unmittelbar benachbarten Teilgebieten erfordern, besitzt diese Methodik eine

ausgesprochen hohe parallele Effizienz (d. h. Beschleunigung gegenüber einer

seriellen Ausführung). Auf den ersten Blick birgt FDS daher ein sehr hohes

Parallelisierungspotenzial in sich.

3.1 DIE FDS-DRUCKGlEICHUNG

Eine Sonderstellung im Parallelisierungsprozess nimmt die Druckgleichung für

den modifizierten Druckterm H ein,

die innerhalb des Gesamtablaufs sehr eng mit der Berechnung aller thermo-

dynamischen Größen verwoben ist. Dabei handelt es sich um eine elliptische

partielle Differentialgleichung vom Poisson-Typ mit verschiedenen Typen von

Randbedingungen, siehe [11].

Die Diskretisierung der Druckgleichung unter Verwendung einer zellzentrierten

Finite-Differenzen-Methode zweiter Ordnung,

führt auf ein globales Gleichungssystem

wobei A eine Matrix des Rnxn bzw. x, b Vektoren des Rn sind. Detaillierte Hinter-

grundinformationen finden sich unter [11, 16].

Bei der Konzeption von lösungsverfahren für die Druckgleichung muss man

sich insbesondere eine sehr spezielle Charakteristik elliptischer partieller

Differentialgleichungen vor Augen führen: Gleichungen dieses Typs besitzen ge-

12


wissermaßen eine unendliche Ausbreitungsgeschwindigkeit für Informationen.

Bereits ein einziger Zeitschritt genügt, um lokale Druckinformationen über das

gesamte Gebiet hinweg zu verteilen, vergleiche Bild 1. lokale Effekte oder Stö-

rungen in den Daten beeinflussen daher unmittelbar die lösung im gesamten

Berechnungsgebiet.

Um eine belastbare parallele Approximation des Druckes zu erzielen, muss auch

der verwendete parallele löser dieses ausgeprägt global orientierte Verhalten

bestmöglich reproduzieren. Im Gegensatz zum sehr lokalen Charakter der

oben genannten expliziten Methodik, erfordert die numerische lösung dieses

Systems daher die Verwendung neuer, gebietsübergreifender Parallelisierungs-

konzepte, die im Mittelpunkt der folgenden Betrachtungen stehen.

3.2 AlGORITHMISCHE BESCHREIBUNG DES AKTUEllEN FFT-DRUCKlÖSERS

Der bisherige Zugang zur lösung der FDS-Druckgleichung basiert auf der Ver-

wendung einer hoch optimierten FFT-Methode aus CRAYFISHPAK [5]. Bedenkt

man, dass FDS in seinen Anfängen nur für rein serielle Single-Mesh-Rechnungen

konzipiert war, so handelte es sich bei dieser Wahl um eine ausgesprochen sinnvolle

Entscheidung, da FFT-Methoden seit vielen Jahren äußerst effizient und erfolgreich

in verschiedensten Wissenschaftszweigen eingesetzt werden. Dies hat sich auch

für Single-Mesh-Rechnungen mit FDS bewahrheitet, für die sich der FFT-Drucklöser

seit Jahren als extrem robust und konkurrenzlos laufzeiteffizient erweist.

Erst die Einführung der Gebietszerlegungstechnik in der Version 4 von FDS ließ

die Verwendung der FFT-Methode im neuen licht erscheinen, da sich immer

wieder Abweichungen zwischen Multi-Mesh-Rechnungen und ihren Single-Mesh-

Entsprechungen feststellen ließen. Dieser Sachverhalt ist wesentlich darauf zu-

Bild 1:

Extrem schnelle

Ausbreitung

der Einström-

Information in der

Druckgleichung

am Beispiel eines

pipeförmigen

Gebietes in 2D

13


rück zu führen, dass FFT-Methoden zur Klasse der direkten Methoden gehören.

Ein weiterer prominenter Vertreter dieser Verfahrensklasse ist beispielsweise

das bekannte Gaußsche Eliminationsverfahren.

Direkte Methoden berechnen die lösung eines Gleichungssystems innerhalb

eines einzigen, möglicherweise sehr komplexen Berechnungszyklus ohne Zwi-

schennäherungen. Die hohe leistungsfähigkeit direkter Methoden basiert dabei

wesentlich auf hoch-rekursiven, gebietsübergreifenden Daten-abhängigkeiten,

die möglichst zeitgleich (häufig in Form kleinster Datenpakete) zur Verfügung

stehen und laufzeiteffizient verarbeitet werden müssen.

Im Falle einer Single-Mesh-Rechnung befinden sich alle hierfür benötigten Daten

im Speicher des einzigen verwendeten Prozessors und können unmittelbar zur

Berechnung heran gezogen werden. Im Falle einer Multi-Mesh-Rechnung werden

die rekursiven Abhängigkeiten jedoch durch die Unterteilung in einzelne Teilgebiete

willkürlich aufgebrochen. Die für die Rechnung erforderlichen Datenpakete sind

über eine Vielzahl von Prozessoren verteilt und müssten zur exakten lösung des

globalen Systems ständig über das gesamte Gebiet hinweg ausgetauscht werden,

was die parallele Effizienz drastisch reduzieren würde. Erschwerend kommt hinzu,

dass bereits vorhandene FFT-Parallelisierungsstrategien in der Regel nicht mit geo-

metrisch motivierten Gebietszerlegungen wie bei FDS konform gehen, siehe [9].

Bedauerlicherweise ist ein Aufbrechen der rekursiven Datenabhängigkeiten zur

Erhöhung des Parallelitätsgrades in der Regel unvermeidbar. Diesem Weg sind

auch die Hauptentwickler gefolgt und haben eine nahe liegende Parallelisie-

rungsstrategie angewandt, die sich einfach in die vorhandene Programmstruk-

tur integrieren ließ: Die Finite-Differenzen-Diskretisierung findet nun nicht mehr

auf globaler Ebene, sondern lokal auf jedem einzelnen Teilgebiet statt, so dass

letztlich n lokale Gleichungssysteme entstehen,

Anstelle der übergreifenden lösung des Drucksystems durch ein einziges,

globales FFT-Verfahren wie im Single-Mesh-Fall, führt nun jedes Teilgebiet ein

eigenes, lokales FFT-Verfahren durch, vergleiche Bild 2. Zur wechselseitigen

Kopplung der lokalen Teilgebietslösungen wird innerhalb jedes Zeitschritts zu-

sätzlich ein koordinierter Datenaustausch zwischen unmittelbar benachbarten

Teilgebieten vorgenommen. Insgesamt reiht sich dieser Ansatz konzeptionell in

die Reihe der additiven Schwarz‘schen Gebietszerlegungsmethoden ein [10].

14


Bild 2:

Parallelisierung

durch verteilte

FFT-Verfahren

mit lokalem

Datenaustausch

im aktuellen FDS-

Drucklöser

Da die Druckgleichung zur Klasse der elliptischen Randwertprobleme gehört,

deren korrekte numerische lösung das Vorhandensein wohldefinierter Rand-

werte voraussetzt, ist man bei diesem Ansatz dazu gezwungen, entlang der neu

entstandenen, inneren Teilgebietsränder möglichst sinnvolle Randbedingungen

vorzugeben. Tatsächlich sind diese jedoch vollkommen unbekannt bzw. müssen

aus bereits berechneten Daten aus dem Gebietsinneren künstlich herleitet

werden.

Mit Blick auf die immanent hohe Ausbreitungsgeschwindigkeit des Drucks wird

schnell klar, dass das unphysikalische Aufbrechen des Gesamtzusammenhangs

durch Verwendung rein lokaler Datenaustausch-Techniken nicht hinreichend

behoben werden kann: Ein Informationsaustausch von einem Gebietsende

zum anderen kann nur schrittweise über alle dazwischen liegenden Teilgebiete

hinweg erfolgen.

Dieser Sachverhalt wird anschaulich in Bild 3 demonstriert, dem eine Untertei-

lung des zuvor betrachteten pipeförmigen Gebietes in vier Teilgebiete zugrunde

liegt: Die Einström-Information vom linken Gebietsrand muss sukzessive durch

das Gesamtgebiet transportiert werden und benötigt insgesamt vier Kommuni-

kationszyklen, um schließlich den rechten Gebietsrand zu erreichen.

Dies steht im krassen Widerspruch zu der tatsächlichen Ausbreitungsgeschwin-

digkeit des Drucks und kann sich umso drastischer auswirken, je höher die An-

zahl der verwendeten Teilgebiete bzw. die damit einhergehende Zerstückelung

des physikalischen Zusammenhangs ist.

15


letztlich kann dieser Sachverhalt zu einer deutlichen Verschlechterung der

Konvergenz beitragen, da globale Abhängigkeiten nicht korrekt wiedergegeben

werden. Trotz extrem hoher paralleler Effizienz können daher für bestimmte

Konstellationen massive Einbrüche in der numerischen Effizienz (Approxi-

mationsgüte bzw. Konvergenzgeschwindigkeit) auftreten, die im schlimmsten

Fall zu numerischer Instabilität und Divergenz führen. Das vorrangige Ziel, die

Skalierbarkeit von FDS auf höhere Prozessorzahlen zuverlässig zu gewährlei-

sten, scheint durch die rein lokale Verwendung von FFT-Verfahren folglich nicht

realisierbar zu sein. Die Hauptentwickler sind sich dieses Sachverhalts bewusst

und bezeichnen die Notwendigkeit zur Verbesserung des Drucklösers in ihrer

Roadmap [2] als vordringliches Ziel:

„We still need to improve the coupling of the pressure solver across mesh

boundaries in a multi-mesh (often called „multi-block“) simulation. Discontinu-

ities in the pressure solution across these boundaries is currently the biggest

speed bump to broader, more reliable use of the multi-mesh approach ... Poten-

tial Research Topic: Improve efficiency of parallel version of FDS on very large

computing clusters (100+ processors) ... Priority: High.“

Bild 3:

Stockender

Datenfluss

vom linken

zum rechten

Gebietsrand beim

parallelen FFT-

Drucklöser

16


Zur lösung dieses Dilemmas ist es unabdingbar, dass das parallele Verfahren

einen Mechanismus zur adäquaten Wiedergabe globaler Effekte beinhaltet.

Eine aktuelle Weiterentwicklung des parallelen FFT-Verfahrens am NIST, das

FFT-Iterationsschema, beruht entsprechend darauf, durch Hinzunahme eines

übergeordneten Grobgitterproblems den gebietsübergreifenden Datenfluss

besser aufzulösen, siehe Bild 4.

Bild 4:

Iterative

Durchführung

der lokalen

FFT-Verfahren in

Kombination mit

einem Grobgitter-

problem im

neuen FFT-

Iterationsschema

Zu diesem Zweck werden die herkömmlichen lokalen FFT-Verfahren in eine

umgreifende Iteration eingebettet und mit einem globalen Grobgitterproblem

gekoppelt. Als Abbruchkriterium für diese Iteration kann im FDS-Eingabefile die

Genauigkeit vorgegeben werden, mit der die Normalkomponenten der lokalen

Geschwindigkeiten entlang innerer Teilgebietsränder übereinstimmen müssen.

Für eine gewünschte Genauigkeit von 10-2 wird z. B. folgender Wert spezifiziert:

&PRES VElOCITY_TOlERANCE = 0.01.

Wie wir später noch sehen werden, lässt sich hiermit tatsächlich eine verbes-

serte Skalierbarkeit für viele Fälle erreichen. Das Konvergenzverhalten dieser

Variante zeigt sich zur Zeit jedoch noch als etwas unzuverlässig bzw. scheint

von der Wahl optimaler Verfahrensparameter abzuhängen, die a-priori nicht

bekannt sind.

17


3.3 AlGORITHMISCHE BESCHREIBUNG DES AlTERNATIVEN SCARC-DRUCK-

lÖSERS

Der alternative Drucklöser ScaRC versucht diese lücke zu schließen und eine

maßgebliche Verbesserung bei der Qualität der Drucklösung zu erzielen. Durch

die ausgeklügelte Kombination von hochauflösender Feingitter-Information mit

gebietsübergreifender Grobgitter-Information auf Basis von iterativen Mehr-

gitter-Techniken wird eine deutlich verbesserte globale Kopplung und Skalier-

barkeit auf große Teilgebietszahlen erreicht. Auf die mathematisch ausführliche

Beschreibung von ScaRC bzw. Mehrgitterverfahren soll an dieser Stelle verzich-

tet werden. Eine detaillierte Herleitung findet sich beispielsweise in [7, 10, 13].

Stattdessen wollen wir uns nachfolgend auf eine informelle Beschreibung der

wichtigsten Eigenschaften dieser Verfahrensklasse beschränken.

Mehrgitterverfahren zählen zu den effizientesten Verfahren für die lösung

großer Gleichungssysteme, mit deren Hilfe sich exzellente, von der Gitterweite

unabhängige Konvergenzraten bei gleichzeitig moderater rechnerischer Kom-

plexität erzielen lassen. Man unterscheidet hierbei zwei verschiedene Zugänge:

• Geometrische Mehrgitterverfahren (GMG) verwenden die Geometrie-Infor-

mation des vorliegenden Problems, um die einzelnen Verfahrensbestand-

teile zu definieren und setzen daher bestimmte geometrische Grundstruk-

turen voraus.

• Algebraische Mehrgitterverfahren (AMG) verwenden lediglich die Informa-

tion, die sich aus dem linearen Gleichungssystem ergibt und sind daher

auch für Probleme auf komplizierteren Berechnungsgebieten bzw. unstruk-

turierten Gittern geeignet.

Beide Varianten sind Bestandteil vieler international verwendeter CFD-Pakete.

Insbesondere AMG-Verfahren finden aufgrund ihrer geometrischen Vielseitig-

keit zunehmend Verwendung, siehe zum Beispiel [8, 12].

Ein wichtiger Unterschied zur Klasse der direkten Verfahren, wie beispielsweise

dem FFT-Verfahren, besteht darin, dass Mehrgitterverfahren zur lösung eines

Gleichungssystems mehrere Berechnungszyklen benötigen. Damit gehören

Mehrgitterverfahren zur Klasse der iterativen Verfahren. Die einzelnen Berech-

nungszyklen sind jedoch von deutlich geringerer Komplexität als dies beim

einzigen Berechnungszyklus eines direkten Verfahrens der Fall ist.

18


Innerhalb dieser Zyklen nähern sich die erzeugten Approximationen immer

mehr der gesuchten lösung an. Insgesamt sind iterative Verfahren häufig ein-

facher zu implementieren als direkte Verfahren, da sie sich im Wesentlichen auf

eine Reihe von Kernkomponenten (Matrix-Vektor-Produkte, Vektor-linear-kom-

binationen und Skalarprodukte) zurückführen lassen, die in geeigneter Weise

optimiert und an die vorliegende Rechnerarchitektur angepasst werden können.

Mehrgitter beruhen im Kern auf der Verwendung eines simplen Defektkorrektur-

Verfahrens

mit einer sogenannten Vorkonditionierungsmatrix B, die grundlegende Eigen-

schaften der Matrix A in sich vereinen sollte, B ∼ part (A), und zur Beschleuni-

gung des lösungsprozesses dient.

Prominente Vertreter dieser Verfahrensklasse sind beispielsweise das Jacobi-

oder Gauß-Seidel-Verfahren, siehe [10, 17]. Gemessen in einer geeigneten Norm

beschreibt der zugehörige Defekt

wie gut die aktuelle Approximation xm das Gleichungssystem Ax = b bereits erfüllt

und stellt somit ein Maß für die Konvergenzgeschwindigkeit des Verfahrens dar.

Eine herkömmliche Basisiteration konvergiert typischerweise jedoch relativ

langsam. Das wesentliche Konstruktionsmerkmal von Mehrgitterverfahren

besteht nun darin, die Konvergenzgeschwindigkeit der Basisiteration zu be-

schleunigen, indem nicht nur auf einem einzigen (feinen) Gitter gearbeitet wird,

sondern eine ganze Hierarchie gröber werdender Gitter zum Einsatz kommt.

Diese Vorgehensweise wird in Bild 5 grafisch veranschaulicht: Ausgehend von

einer Startlösung auf dem feinsten Gitterlevel werden zunächst einige Schritte der

Basisiteration durchgeführt. Auf diese Weise wird in der Regel bereits nach we-

nigen Iterationen eine deutliche Reduktion der hochfrequenten Fehleranteile des

Defektes erreicht. Diese besondere Eigenschaft der Basisiteration wird mit dem

Begriff Glättungseigenschaft bezeichnet. Die niederfrequenten Fehleranteile blei-

ben jedoch üblicherweise nahezu unverändert und können noch sehr groß sein.

Dieser Sachverhalt legt es nahe, den so geglätteten Defekt auf ein vergröbertes

Gitter zu restringieren, d. h. durch geeignete Mittelwerte abzubilden. Die nieder-

19


Bild 5:

Mehrgitter-

Hierarchie mit

verschieden

groben Gittern

sowie Restriktion

und Prolongation

des Defektes

frequenten Anteile des feinen levels erscheinen auf dem gröberen level wiede-

rum als hochfrequente Anteile und können mit Hilfe einer weiteren Basisiteration

erneut geglättet werden, hier jedoch zu deutlich geringeren Berechnungskosten.

Diese Vorglättung wird rekursiv so lange fortgesetzt, bis ein gröbster Gitterlevel

erreicht ist, auf dem das verbleibende Grobgitterproblem mit sehr geringem

Aufwand exakt gelöst werden kann. Anschließend wird die erzielte Grobgitterlö-

sung mit Hilfe geeigneter Interpolationsmechanismen sukzessive auf die feineren

Gitter zurück prolongiert, wobei wiederum auf jedem level einige Schritte der

Basisiteration zur Nachglättung durchgeführt werden können.

Jeder Gitterlevel ist für die Reduktion eines ganz bestimmten Frequenzbe-

reiches der Fehleranteile zuständig. Die Effizienz der gesamten Methode basiert

wesentlich darauf, wie gut diese Bereiche interagieren bzw. die verschiedenen

Gitterweiten aufeinander abgestimmt sind. Insgesamt tragen das Durchlaufen

der Gitterhierarchie sowie die lösung des Grobgitterproblems zu einer ausge-

sprochen starken gebietsübergreifenden Kopplung bei.

Wie aus Bild 5 hervor geht, ist auf jedem Gitterlevel nur lokale Kommunikation

zwischen benachbarten Teilgebieten für die Berechnung der benötigten Matrix-

Vektor-Produkte erforderlich. lediglich die Berechnung des Grobgitterproblems

erfordert die weniger effiziente Durchführung globaler Kommunikation.

Zur Wahl der Vorkonditionierungsmatrix können nun auf sehr natürliche Weise

Gebietszerlegungsstrategien herangezogen werden, beispielsweise

20


Die Vorkonditionierung besteht folglich aus der lösung der einzelnen Teilgebiets-

probleme, wodurch sehr viele spezifische Eigenschaften des Gesamtsystems

ausgenutzt werden können.

Aufbauend auf diesen iterativen Kernbestandteilen definiert sich der verallge-

meinerte Gebietszerlegungslöser ScaRC nun als fein abgestimmte Symbiose aus

effizienten globalen und lokalen Techniken, die in unterschiedlichen Varianten dem

vorliegenden Problem bestmöglich angepasst werden können. Hierzu versucht

ScaRC in geeigneter Weise gebietsübergreifende Mehrgitteransätze mit hoch-

effizienten lokalen lösern für die Vorkonditionierungsprobleme innerhalb der

einzelnen Teilgebiete zu kombinieren. Dabei handelt es sich nicht um ein einziges

fest gefügtes Verfahren, sondern vielmehr um eine Art Baukastensystem, das ein

großes Spektrum moderner Gebietszerlegungs- und Mehrgitterzugänge für die lö-

sung partieller Differentialgleichungen enthält und sinnvoll kombiniert. Detaillierte

Information zur algorithmischen Herleitung von ScaRC finden sich in [13].

Ein fundamentales Konstruktionsmerkmal von ScaRC besteht darin, ein möglichst

hohes Maß an Datenlokalität zu bewahren, um lokale Prozessorbestandteile wie

z. B. den Cache bestmöglich ausnutzen zu können. lokale Irregularitäten können

innerhalb einzelner Teilgebiete „versteckt“ und dort mit besonders robusten

lokalen lösern oder feineren Gittern behandelt werden. Wie in Bild 6 illustriert,

könnte dies zukünftig den Weg ebnen hin zur adaptiven Verfeinerung lokaler Be-

reiche mit besonders hoher Variationen der betrachteten physikalischen Größen

auf Basis von Techniken des Adaptive Mesh Refinement (AMR).

Der größte Vorteil dieses Zugangs liegt darin, dass ScaRC im Gegensatz zum

parallelen FFT-Ansatz von einer globalen Diskretisierung des Gesamtgebietes

ausgeht. Die zugehörige globale Matrix A innerhalb der Basisiteration ist jedoch

ein rein formaler Ausgangspunkt und wird nie zusammenhängend aufgebaut,

sondern über die einzelnen Prozessoren verteilt gespeichert. Die gesamte

Iteration ist global definiert und verwendet die lösung der Teilgebietsprobleme

lediglich zur näherungsweisen Korrektur des globalen Defektes. Alle globalen

Matrix-Vektor- und Vektor-Operationen liefern dasselbe Ergebnis wie bei einer

hypothetischen seriellen Berechnung des Gesamtsystems, was als datenparallele

Ausführung bezeichnet wird.

Dies bedeutet insbesondere, dass keine künstlichen Bedingungen entlang innerer

Gebietsränder gestellt werden müssen, da innere Randpunkte - bezogen auf die

virtuelle globale Matrix A - als ganz normale innere Punkte behandelt werden.

21


Bild 6:

lokale adaptive

Verfeinerung

einzelner

Bereiche des

Gebietes

Am Ende der Berechnung steht eine globale lösung, die entlang innerer Gebiets-

kanten vollkommen konsistent ist, ganz im Gegensatz zu der eher lockeren Zusam-

mensetzung der lokalen lösungen innerhalb des herkömmlichen FFT-Zugangs bzw.

zum neuen FFT-Iterationsschema, das nach Konstruktion nur eine vorgegebene

Genauigkeit für die Geschwindigkeitskomponenten entlang innerer Ränder erzielt.

3.4 ENTWICKlUNGSSTAND DES SCARC-DRUCKlÖSERS

Wie bereits angedeutet, handelt es sich bei ScaRC nicht um ein einziges

Verfahren, sondern um eine umfangreiche Verfahrensklasse, die neben den

genannten Mehrgittertechniken noch weitere bewährte iterative lösungstech-

niken beinhaltet. So steht auch ein datenparallelisiertes Konjugierte Gradienten

Verfahren (CG) bzw. seine Erweiterung für nicht-symmetrische lineare Systeme,

das stabilisierte Bikonjugierte Gradienten Verfahren (BICGstab) zur Verfügung.

Details hierzu finden sich beispielsweise unter [17].

Alle nachfolgend präsentierten Testrechnungen basieren auf der Verwendung eines

geometrischen Mehrgitterverfahrens mit SSOR-Glättung, siehe [10, 17]. Die geome-

trische Mehrgittervariante zeichnet sich durch eine besonders hohe numerische

Effizienz mit ausgesprochen niedrigen Konvergenzraten aus. Da die gröberen Gitter

jeweils durch Verdopplung der Gitterweite des nächst feineren Gitters erzeugt wer-

den, entstehen auf den einzelnen Prozessoren sehr reguläre Datenstrukturen, zu de-

ren Behandlung hoch-optimierte Bibliotheken herangezogen bzw. lokale Prozessor-

Features wie beispielsweise der Cache optimal ausgenutzt werden können.

leider kann diese „geometrische“ Konstruktionsweise aber auch eine Ein-

schränkung bzgl. der potentiell verwendbaren Gitterweiten darstellen, da die

Zellanzahlen in den einzelnen Raumdimensionen notwendigerweise Zweier-

potenzen sein müssen. Dieser Umstand kann jedoch dadurch abgemildert

22


werden, dass man Zellanzahlen verwendet, die zwar keine Zweierpotenz sind,

aber zumindest einige Male durch zwei geteilt werden können. Dies reduziert

letztlich die Anzahl der resultierenden Gitterlevel. Unter Umständen wird nicht

das gröbst mögliche Grobgitter (nämlich die Gebietszerlegung selbst) erreicht,

so dass die exakte lösung des Grobgitterproblems sich aufwändiger gestalten

kann. Da die Gitterlevel lediglich zur näherungsweisen Korrektur des Defekts

eingesetzt werden, kann aber auch im Fall einer unvollständigen Gitterhierar-

chie in der Regel eine sehr gute Approximation erreicht werden.

Um jedoch völlige Freiheit in der Wahl der möglichen Gitterweiten und Zell-

anzahlen zu erreichen, arbeiten wir zur Zeit intensiv an der Erweiterung von

ScaRC auf die Klasse der algebraischen Mehrgitterverfahren. Zur Erzeugung der

gröberen Gitter werden dabei lediglich diejenigen Informationen verwendet die

in der Matrix selbst vorhanden sind. Es wird jedoch keine feste Gitterstruktur

vorausgesetzt. Im Gegenteil: Zur Erzeugung der gröberen Gitter stehen eine

ganze Reihe verschiedenartiger Vergröberungstechniken zur Verfügung, die auf

sehr unterschiedliche „grobe“ Gitterstrukturen führen können, vergleiche [8].

Insbesondere kann auf den gröberen Gittern unter Umständen keine zeilenwei-

se Nummerierung der Gitterpunkte mehr vorgenommen werden.

Für Single-Mesh-Rechnungen wurde die AMG-Methodik bereits erfolgreich

implementiert. Aufgrund der Vielzahl an möglichen Strukturen für die gröberen

Gitter, muss für Multi-Mesh-Fälle eine umfassende Erweiterung der bereits vor-

handenen Techniken zum Datenaustausch zwischen benachbarten Teilgebieten

vorgenommen werden.

Gleichzeitig können die resultieren Gitter bzw. Matrizen nicht mehr mit den in

FDS vorhandenen Datenstrukturen erfasst werden, da diese speziell auf zeilen-

weise nummerierte Gitter abgestimmt sind. Daher wurde in jüngster Vergangen-

heit eine neue Technik zur Speicherung der entstehenden Gitter- und Matrix-In-

formationen in FDS implementiert. Diese sogenannte kompakte Speichertechnik

birgt ein hohes Maß an Variabilität in sich. Ein besonderer Vorteil dieser Technik

besteht insbesondere darin, dass sie für die Speicherung adaptiv verfeinerter

Gitter sehr gut geeignet ist.

Die Arbeiten am parallelen AMG sind zum gegenwärtigen Zeitpunkt bereits weit

fortgeschritten. Wir hoffen sehr, damit mittelfristig eine allgemeingültige ScaRC-

Version zur Verfügung stellen zu können, die zur nachhaltigen Verbesserung

der FDS-Drucklösung beitragen wird.

23


Die Ergebnisse erster numerischer Testreihen auf Basis von GMG-Verfahren stim-

men sehr hoffnungsvoll: Wie die im Anschluss dargestellten Resultate belegen, wird

eine hohe physikalische Genauigkeit und eine deutlich verbesserte Skalierbarkeit

im Vergleich zum aktuellem FFT-Drucklöser erzielt. In Absprache mit den Hauptent-

wicklern wurde daher bereits eine Beta-Version von ScaRC mit speziell angepassten

Schnittstellen in das offizielle FDS-Repository integriert. Die Implementierung der

zugehörigen Komponenten baut wesentlich auf bereits vorhandenen Datenstruk-

turen auf, sodass die herkömmliche Eingabelogik, die weltweit einer Vielzahl von

Benutzern vertraut ist, erhalten bleibt. In der kommenden FDS Version 6 soll ScaRC

dann als optionaler Drucklöser offiziell zur Verfügung stehen.

4 NUMERISCHE TESTRECHNUNGEN

Von einem korrekten parallelen lösungsalgorithmus wird erwartet, dass er

unabhängig von der Anzahl der Teilgebiete immer die gleiche oder doch nahezu

die gleiche lösung produziert. Wie bereits erwähnt, kann es vor dem Hinter-

grund der notwendigen Modifikationen zur Erhöhung der parallelen Effizienz,

die sich aus dem Aufbrechen rekursiver Strukturen ergeben, unter Umständen

erforderlich sein, leichte Abweichungen zu akzeptieren. Diese müssen jedoch in

jedem Fall möglichst klein und kontrollierbar bleiben.

Nachfolgend soll am Beispiel zweier Testreihen untersucht werden, inwieweit

die betrachteten Drucklöser dieser Anforderung genügen. Zur Vereinfachung

der Notation werden dabei folgende Bezeichnungen verwendet:

• FFT: altes FFT-Verfahren mit rein lokaler lösung der Teilgebietsprobleme,

• FFTI, tol=10-n: neues FFT-Iterationsschema mit lokaler lösung der Teil-

gebietsprobleme plus zusätzlichem Grobgitterproblem im Rahmen einer

übergreifenden Iteration für eine vorgegebene Genauigkeitstoleranz von

10-n bezüglich der Geschwindigkeitskomponenten entlang innerer Ränder,

• ScaRC: neues ScaRC-Verfahren mit lösung der Teilgebietsprobleme auf

einer Gitterhierarchie plus Grobgitterproblem.

Bei der Konzeption der Testreihen wurde großen Wert darauf gelegt, die Effekte

des Drucklösers so gut wie möglich von den Effekten der anderen lösungs-

komponenten zu separieren, um letztlich aussagekräftige Rückschlüsse auf die

Qualität des reinen Drucklösers ziehen zu können.

24


Aufgrund sehr komplexer Wechselwirkungen der einzelnen Verfahrensbestand-

teile ist es in der Regel ausgesprochen schwierig, die Auswirkungen einzelner

Komponenten für sich alleine zu betrachten. Im ungünstigsten Fall kann es zu

schwer nachvollziehbaren Überlagerungen oder gegenseitigen Auslöschungen

von Fehlereinflüssen kommen, die eine objektive Bewertung erschweren. Detail-

lierte Konzepte für ein umfassendes Verifikationskonzept zur Eignungsprüfung

von FDS wurden in Kooperation mit Matthias Münch, INURI GmbH, erarbeitet

und sind Grundlage der folgenden Testreihen, siehe [13, 14].

Alle nachfolgend dargestellten Testrechnungen wurden mit einer aktuellen

Entwicklerversion, FDS 5.5.3 Revision 8305, durchgeführt und spiegeln somit

den gegenwärtigen Entwicklungsstand von FDS wider. Zur Kontrolle wurden die

dargestellten Rechnungen auch mit der offiziell verfügbaren Programmversion,

FDS 5.5.3 Revision 7031, verglichen. Dabei konnten aber keine nennenswerten

qualitativen Unterschiede bezüglich der Drucklösung festgestellt werden.

4.1 DURCHSTRÖMUNG EINES KANAlS IN 2D

Die folgende Testreihe basiert auf einer einfachen Kanal-Geometrie in 2D. Von

der linken Seite findet eine Einströmung mit vorgegebener Geschwindigkeit

statt, die rechte Seite des Gebietes ist offen. Um die Skalierbarkeit der betrach-

teten Drucklöser miteinander zu vergleichen, wurde das Gebiet wie in Bild 7

skizziert in N = 4, 8 und 16 Teilgebiete in vertikaler Richtung unterteilt.

Die zugehörige analytische lösung lautet in diesem Fall

Eine detaillierte Beschreibung dieser Testreihe findet sich unter [14]. Im Fall

einer linearen Einströmung,

ergibt sich für die Druckdifferenz zwischen Ein- und Ausströmrand eine analytische

Druckdifferenz von

was einen unmittelbaren Vergleich mit den numerisch ermittelten Daten erlaubt.

Den Testrechnungen liegen folgende Werte zugrunde: u0=1 m⁄s , ∆x=0.8,

ρ=1.19882 kg⁄m3 sowie ein Endzeitpunkt von T=4 s. Die Eingabeparameter

25


Bild 7:

Durchströmung

eines Kanals in 2D

mit vorgegebener

Einströmung

und zugehöriger

Gebietszerlegung

in vertikale

Streifen

wurden so gewählt, dass weder Strahlungs- noch Reibungseffekte auftraten

und ein konstanter Hintergrunddruck vorlag. Außerdem wurde eine direkte

numerische Simulation (DNS) durchgeführt (der FDS-Turbulenzlöser wurde

ausgeschaltet, alle potentiell auftretenden turbulenten längenskalen wurden

berechnet).

Sehr deutlich wird der ausgesprochen lokale Charakter der herkömmlichen,

parallelen FFT-Methodik zunächst am Beispiel einer Geschwindigkeitsdarstel-

lung im 8-Mesh-Fall, vergleiche Bild 8. Während ScaRC die Einström-Informati-

on unmittelbar durch das gesamte Gebiet transportiert und damit das tatsäch-

liche physikalische Verhalten reproduziert, wird die Information im FFT-Fall nur

stockend übertragen.

Bild 8:

Stockender

Informations-

transfer im

parallelen FFT-

Verfahren für den

8-Mesh-Fall

Weiterhin stellen die Bildersequenzen 9, 10 und 11 für alle betrachteten Teilge-

bietszahlen die jeweils numerisch ermittelten Druckdifferenzen der analytischen

Druckdifferenz gegenüber.

26


Ganz offensichtlich sind alle verschiedenen Drucklöser im Single-Mesh-Fall dazu

in der lage, die analytische Druckdifferenz exakt zu reproduzieren (gerade

linie), was für die Korrektheit des Single-Mesh-Zugangs spricht. Betrachtet man

jedoch das Standard-FFT-Verfahren in Bild 9, so zeigen sich deutliche Abwei-

chungen für die Multi-Mesh-Fälle, die sich umso stärker bemerkbar machen, je

mehr Teilgebiete verwendet werden. Das lokale Vorgehen ist offensichtlich nicht

in der lage, den globalen Druckverlauf adäquat zu erfassen.

Auch im Fall des FFT-Iterationsschemas kann für keine der betrachteten Tole-

ranzen eine völlige Übereinstimmung erreicht werden, siehe Bild 10. Während

die Fälle tol = 10-1 und 10-2 keinen qualitativ nennenswerten Unterschied zum

Standard-FFT-Verfahren bewirken, zeigen die Fälle tol = 10-3 und insbesondere

10-4 tatsächlich Verbesserungen, die jedoch nur von kurzer Dauer sind. Auch

der anfänglich gute Verlauf im Fall tol = 10-4 wird zunehmend oszillativ und

schaukelt sich mit fortschreitender Berechnungsdauer auf. Die Anzahl der

benötigten Iterationen des FFT-Iterationsschemas zum Erreichen der Geschwin-

digkeitstoleranz von 10-4 schwankt dabei ganz erheblich zwischen 7 und 577

Iterationen! Damit einher geht eine Vervierfachung der Rechenzeit im Vergleich

zum herkömmlichen parallelen FFT-Verfahren. Trotzdem ist kein befriedigendes

Konvergenzverhalten erkennbar.

Dagegen reproduziert ScaRC die analytisch ermittelte Druckdifferenz gleicher-

maßen für alle betrachteten Multi-Mesh-Fälle ohne jegliche Abweichungen,

vergleiche Bild 11. Dieses Resultat ist außerdem unabhängig von der Art der be-

trachteten Gebietszerlegung. Zum Vergleich wurden entsprechende Tests auch

für Gebietszerlegungen mit bis zu 64 Teilgebieten nur in horizontaler Richtung

sowie in vertikaler und horizontaler Richtung durchgeführt. Auch für diese Fälle

konnte eine völlige Übereinstimmung erzielt werden.

Insgesamt illustriert dieses Testbeispiel sehr eindrucksvoll die Vorteile eines

Verfahrens, dessen numerische Konzeption den Eigenheiten des zugrunde

liegenden Problems angepasst ist.

27


Bild 9:

Vergleich der

Single- und

Multi-Mesh-

Druckdifferenzen

mit der

analytischen

Druckdifferenz für

das herkömmliche

FFT Verfahren

Bild 11:

Vergleich der

Single- und

Multi-Mesh-

Druckdifferenzen

mit der

analytischen

Druckdifferenz

für das ScaRC-

Verfahren

Bild 10:

Vergleich der

Single- und

Multi-Mesh-

Druckdifferenzen

mit der

analytischen

Druckdifferenz

für das FFT-

Iterationsschema

mit

unterschiedlichen

Toleranzen für die

Geschwindigkeits-

komponenten

28


4.2 DREIFACH-EINSTRÖMUNG IN EIN KUBISCHES GEBIET IN 3D

Der nachfolgende Testfall behandelt ein einfaches würfelförmiges Gebiet in 3D,

in das von drei Seiten eine Einströmung mit einer konstanten Geschwindigkeit

von 2 m/s stattfindet. Bild 12 zeigt die zugrunde liegende Geometrie und eine

Geschwindigkeits-Ansicht im Single-Mesh-Fall.

Um einen Eindruck zu erhalten, wie gut sich die betrachteten Drucklöser auf

große Teilgebietszahlen skalieren lassen, wurde das Gebiet in 64, 128 und 512

Teilgebiete unterteilt, vergleiche Bild 13. Die Verwendung hunderter von Teil-

gebieten – ein erklärtes Ziel der Hauptentwickler [2] – ist in der gegenwärtigen

Anwendungspraxis von FDS noch unüblich und stellt quasi einen Härtetest dar.

Bild 12:

Kubisches

Gebiet in 3D mit

Einströmung von

drei Seiten

Bild 13:

Massiv parallele

Gebietszerle-

gungen für ein

würfelförmiges

3D-Gebiet

Da für diesen Fall keine analytische lösung vorliegt, müssen sich die Multi-

Mesh-Resultate aller betrachteten Drucklöser nachfolgend mit dem zugehörigen

Single-Mesh-Resultat messen. In der Bildersequenz 14, 15 und 16 stellt der Ver-

29


Bild 14:

Vergleich der

Multi-Mesh-

Druckverläufe für

das herkömmliche

FFT-Verfahren

Bild 15:

Vergleich der

Multi-Mesh-

Druckverläufe

für das FFT-

Iterationsschema

mit unterschied-

lichen Toleranzen

für die

Geschwindigkeits-

komponenten

Bild 16:

Vergleich der

Single- und

Multi-Mesh-

Druckverläufe

für das ScaRC-

Verfahren

30


lauf des numerisch berechneten Single-Mesh-Drucks im Mittelpunkt des Gebiets

eine Referenzkurve dar, mit der die Verläufe aller Multi-Mesh-Rechnungen

verglichen werden.

Wie aus Bild 14 hervorgeht, zeigt das herkömmliche FFT-Verfahren bereits im

64-Mesh-Fall keine Übereinstimmung mehr mit dem 1-Mesh-Fall. Noch ausge-

prägter sind die Unterschiede für 128 bzw. 512 Teilgebiete, da das rein lokale

Vorgehen der hohen Anzahl an Teilgebieten nicht gerecht wird. Offenbar wer-

den die Druckverläufe der Multi-Mesh-Fälle grundsätzlich eher „unterschätzt“.

Im Fall des FFT-Iterationsschemas in Bild 15 wird für die Toleranz 10-1 noch kein

nennenswerter Unterschied zum herkömmlichen FFT-Verfahren erreicht. Dage-

gen führt der Fall 10-2 tatsächlich zu einer verbesserten Konsistenz zwischen

dem Single-Mesh und den Multi-Mesh-Verläufen. Für Toleranzen kleiner als 10-3

ist erstaunlicherweise jedoch wieder eine deutliche Verschlechterung erkenn-

bar. Insgesamt neigt das FFT-Iterationsschema offenbar dazu, die Verläufe zu

„überschätzen“.

Die besten Resultate können auch hier wieder mit Hilfe des ScaRC-Verfahrens

erreicht werden. In Bild 16 ist eine hohe Übereinstimmung aller betrachteten

Multi-Mesh-Verläufe mit dem Single-Mesh-Verlauf erkennbar. Mit fortschreiten-

der Berechnungsdauer wird jedoch offensichtlich, dass auch für ScaRC zurzeit

noch keine vollständige Konsistenz vorliegt. Je nach Problemkonstellation

können über die Gesamtrechenzeit hinweg immer noch mehr oder weniger

starke Abweichungen entstehen. Dies haben weitere Testrechnungen mit unter-

schiedlichen Eingabeparametern belegt. Auf diesen Sachverhalt werden wir an

späterer Stelle noch einmal zurückkommen.

Im Vergleich zum vorangehenden Testfall zur Durchströmung eines pipeförmigen

Kanals bringt das sonderbare Verhalten des FFT-Iterationsschemas leider einige

Verwirrung mit sich: Während für die Kanaldurchströmung eine Toleranz von

10-2 offensichtlich noch keine hinreichende Verbesserung bewirkt hat und erst

für 10-4 eine deutliche Annäherung zu erkennen war, ist das Verhalten hier ge-

nau umgekehrt. Diese Unsicherheit lässt es sehr schwierig erscheinen, vorweg

einen geeigneten Toleranzwert für eine gegebene Rechnung festzulegen.

Ein näherer Blick auf das Konvergenzverhalten des FFT-Iterationsschemas

bringt hier etwas mehr Klarheit: Im Fall von 128 Teilgebieten mit einer vorgege-

benen Toleranz von 10-3 für die Geschwindigkeitskomponenten stellt Bild 17 die

tatsächlich pro FDS-Zeititeration erreichte Fehlertoleranz der dazu benötigten

31


Anzahl von FFT-Iterationen gegenüber. Es geht deutlich hervor, dass die ange-

strebte Toleranz nur in weniger als der Hälfte aller Zeititerationen tatsächlich

erreicht wird! Mit zunehmendem Rechnungsverlauf wird immer häufiger die

Maximalanzahl von 100 FFT-Iterationen durchgeführt. Dennoch wird gerade in

diesen Fällen in der Regel nur eine Genauigkeit von durchschnittlich 5•10-2 für

die inneren Geschwindigkeitskomponenten erreicht.

Bild 17:

Gegenüberstellung

der erzielten

Genauigkeit für die

Geschwindigkeits-

komponenten

mit der dazu

benötigten Anzahl

an FFT-Iterationen

pro FDS-

Zeititeration

Die Maximalanzahl möglicher FFT-Iterationen wurde für alle hier gezeigten

Rechnungen auf 100 gesetzt. Obwohl die entsprechende Standard-Einstellung,

MAX_PRESSURE_ITERATIONS im &PRES-Block, aktuell bei 10.000 liegt, hat es

sich im Verlauf vieler Testreihen heraus gestellt, dass ein derart hoher Wert kei-

nerlei Verbesserung mit sich bringt, sondern nur zu einem drastischen Anstieg

der Rechenzeit führt. Optimale Einstellungen für die Maximalanzahl der FFT-

Iterationen hängen letztlich von der Gitterweite ab und sind ebenfalls a-priori

nur schwer bestimmbar.

32


Betrachten wir für die letzte durchgeführte Zeititeration im obigen 128-Mesh-

Fall, wie sich die Geschwindigkeitskomponenten entlang innerer Gebietskanten

im Verlauf des FFT-Iterationsschemas einander angleichen, so wird deutlich,

dass relativ schnell eine Reduktion des Fehlers bis auf etwa 6•10-3 erreicht wird,

siehe Bild 18. Danach stagniert das Verfahren jedoch und ist nicht in der lage,

die angestrebte Toleranz von 10-3 zu erreichen. In diesem Fall hätte auch die

sehr viel kleinere Maximalanzahl von 10 FFT-Iterationen genügt, um letztlich

zum selben finalen Ergebnis der Gesamtiteration zu gelangen.

Bild 18: Typische

Fehlerreduktion für die

Geschwindigkeitskomponenten im

FFT-Iterationsschema

Offenbar kann die tatsächlich erzielbare Toleranz für die Geschwindigkeitskom-

ponenten entlang innerer Gebietskanten nicht vollständig vom Anwender selbst

bestimmt werden, sondern hängt von bestimmten Eigenschaften des Problems

bzw. den vorliegenden Diskretisierungsparametern ab. Dieser Zusammenhang

soll zukünftig noch näher analysiert werden.

4.3 DISKUSSION DER TESTREIHEN

Die dargestellten Testreihen belegen eindeutig den Zugewinn an Genauigkeit

und Skalierbarkeit durch die neue ScaRC-Methodik. Da ScaRC im Gegensatz zum

herkömmlichen FFT-Ansatz auf der lösung des globalen Gleichungssystems ba-

siert, werden grundlegende Problemeigenschaften wie insbesondere die starke

globale Kopplung des Drucks gut reproduziert.

leider benötigt der aktuelle ScaRC-löser im Vergleich zum herkömmlichen FFT-

löser noch eine höhere Berechnungszeit (etwa um den Faktor 2-3). Im Fall von

33


Single-Mesh-Rechnungen, wo beide Zugänge auf exakt dieselbe lösung führen,

sollte demnach grundsätzlich das FFT-Verfahren eingesetzt werden. Es muss

allerdings betont werden, dass bisher noch keine größeren Anstrengungen

zur laufzeit-Optimierung von ScaRC vorgenommen wurden, da zunächst die

Ausarbeitung eines stabilen numerischen Grundkonzeptes im Mittelpunkt der

Bemühungen stand. Ausführliche laufzeituntersuchungen stehen daher noch

aus, werden aber zu gegebener Zeit folgen.

Der Vorteil der ScaRC-Methodik macht sich grundsätzlich erst im Multi-Mesh-

Fall bemerkbar. Zwar wird mit dem parallelen FFT-löser schneller der festge-

legte Endzeitpunkt der Berechnungen erreicht (sofern überhaupt Konvergenz

vorliegt). Jedoch unterschieden sich die erzielten Resultate unter Umständen

erheblich von der entsprechenden Single-Mesh-lösung und sind daher von

fragwürdiger Bedeutung. Der Blick auf die reine Rechenzeit kann folglich kein

alleiniges Bewertungskriterium für die Qualität des parallelen lösers sein.

Dieser muss sich vielmehr in erster linie an seinen numerischen Konvergenzei-

genschaften messen lassen. letztlich muss die numerische Effizienz des lösers

im Zweifelsfall höher bewertet werden als die parallele Effizienz.

Um tatsächlich eine konsistente lösung zu erhalten, muss das Multi-Mesh-Ver-

fahren notwendigerweise dasselbe Gleichungssystem lösen wie das Single-Mesh-

Verfahren. Das herkömmliche parallele FFT-Verfahren koppelt jedoch lediglich

die lösung lokaler Gleichungssysteme und wird diesem Anspruch nicht gerecht.

Durch die Verwendung eines Grobgitterproblems kommt das FFT-Iterations-

schema der lösung des globalen Systems schon deutlich näher, jedoch nur bis

zu einer bestimmten Toleranz, die sich nicht beliebig verkleinern lässt. Dagegen

basiert ScaRC nach Konstruktion auf einer globalen Vorgehensweise und ist

dazu in der lage, volle Konsistenz zu erreichen.

Die Problematik, dass auch für ScaRC-Testläufe noch Inkonsistenzen festgestellt

werden können, wurde bereits ausführlich mit den Hauptentwicklern diskutiert

und unterliegt einer dauerhaften Beobachtung. Um ausschließen zu können,

dass die vorliegenden Inkonsistenzen eventuell doch auf Fehler im ScaRC-Code

zurückzuführen sind, wurde anhand verschiedener Testszenarien jeweils der fol-

gende Versuch vorgenommen: Für eine gegebene Konstellation wurde zunächst

der zugehörige Single-Mesh-Fall berechnet. Die dabei ins Gleichungssystem

eingehenden Daten wurden dann entsprechend der vorliegenden Teilgebiets-

struktur „von Hand“ aufgesplittet und der zugehörigen Multi-Mesh-Rechnung

eingespeist. Auf diese Weise war sicher gestellt, dass sowohl die Single- als auch

34


die Multi-Mesh-Rechnung jeweils von exakt denselben Eingabedaten ausgehen.

Und tatsächlich konnte unter diesen Umständen mit ScaRC jeweils völlige Konsi-

stenz zwischen beiden Fällen erzielt werden, was auf ein korrektes datenparalleles

Arbeiten von ScaRC schließen lässt.

Es hat sich inzwischen herausgestellt, dass die beobachteten Inkonsistenzen auf

den „umliegenden“ Code zurück zu führen sind, d. h. auf die Programmbestand-

teile, die den Druck-Code umfassen und ihn mit Eingabedaten wie den äußeren

Randbedingungen und der rechten Seite des Gleichungssystems versorgen. Die-

ser Sachverhalt betrifft alle betrachteten Drucklöser gleichermaßen. In Zusam-

menarbeit mit Kevin McGrattan konnten bereits einige Problemstellen behoben

werden, offensichtlich liegen aber noch weitere Passagen im Code vor, die zu

einer unterschiedlichen Behandlung von Single- und Multi-Mesh-Fällen führen.

Die Verwendung der ScaRC-Methodik macht es erstmalig möglich, die erwähnten

Inkonsistenzen zwischen Single- und Multi-Mesh-Rechnungen überhaupt aufzu-

decken. Die bis dato verwendete FFT-Methodik liefert im Multi-Mesh-Fall bereits

nach Konstruktion NICHT dieselbe lösung wie im Single-Mesh-Fall, sondern stellt

lediglich eine (im günstigsten Fall) gute Näherung dazu dar. Dies erschwert es

ungemein, verschiedene Fehlereinflüsse voneinander zu separieren. Ist jedoch,

wie im Fall von ScaRC, sichergestellt, dass gleiche Eingabewerte im Single- und

Multi-Mesh-Fall tatsächlich auch zur gleichen lösung führen, so lassen sich Fehler-

einflüsse anderer Komponenten einfacher identifizieren.

Es ist mit den Hauptentwicklern vereinbart, in Kürze weitere gemeinsame

Verifikationsfälle auszuarbeiten, die mögliche Fehlerquellen aufdecken und

diese Inkonsistenzen restlos ausmerzen sollen. Sicherlich wird auch das aktuelle

FFT-Iterationsschema von dieser Entwicklung profitieren. Dennoch ist nach

Konstruktion nur ScaRC in der lage, unabhängig von der Teilgebietszahl die

lösung des korrekten globalen Systems zu reproduzieren. Der aktuelle Stand

der weiteren Entwicklungen kann fortlaufend unter [3] eingesehen werden.

35


5 FAZIT

Abschließend ist festzustellen, dass es sich bei FDS um ein ausgesprochen

leistungsfähiges Programm zur numerischen Brand- und Rauchgas-Simulation

handelt, das seine große leistungsfähigkeit bisher insbesondere für Single-

Mesh-Rechnungen entfaltet. Vor dem Hintergrund der wachsenden Bedeutung

des Parallel Computing sind die aktuell verwendeten Parallelisierungskonzepte

ein Schritt in die richtige Richtung, weisen jedoch noch ein deutliches Verbesse-

rungspotential auf. In Kooperation mit den Hauptentwicklern werden von

hhpberlin daher große Anstrengungen unternommen, im Rahmen der weiteren

ScaRC-Entwicklung die hohe leistungsfähigkeit der Single-Mesh-Version auch

für den massiv parallelen Multi-Mesh-Fall nutzbar zu machen. Auf dem Weg

dahin zeichnen sich bereits sehr vielversprechende Erfolge ab, die ein weiteres

Vorgehen rechtfertigen. Inzwischen wurde eine solide Basis für einen alter-

nativen Drucklöser geschaffen, der das Potential für weitere Optimierungen,

insbesondere in Hinblick auf adaptive Gitterverfeinerungen, in sich birgt.

37


lIteratur

[1] http://www.fds-smv.net: Fire Dynamics Simulator, Official Website, Hosted

at the National Institute of Standards and Technology.

[2] http://code.google.com/p/fds-smv/wiki/FDS_Road_Map: The FDS/Smoke-

view Research Plan.

[3] http://code.google.com/p/fds-smv/issues/detail?id=1132: Issue Tracker of

FDS-SMV website, Issue 1132, ScaRC development.

[4] http://www.FDS Usergroup.de: Offizielle Website der FDS Usergroup.

[5] http://www.cisl.ucar.edu/css/software/fishpack.html: Efficient Fortran

subprograms for the Solution of Separable Elliptic Partial Differential

Equations.

[6] http://www.mcs.anl.gov/research/projects/mpi: The Message Passing

Interface standard.

[7] Briggs, W. l., Henson V. E., McCormick S. F.: A Multigrid Tutorial, SIAM,

2000.

[8] Bruaset A. M., A.Tveioto A (eds.): Numerical Solution of Partial Differential

Equations on Parallel Computers, lecture Notes in Computational Science

and Engineering, Springer, 51, 2006.

[9] Frommer, lösung linearer Gleichungssysteme auf Parallelrechnern, Fried-

rich Verlag & Sohn, 1990.

[10] Hackbusch, W.: Iterative lösung großer schwachbesetzter Gleichungs-

systeme, Teubner Studienbücher Mathematik, Stuttgart, 1993.

[11] Hackbusch, W.: Theorie und Numerik elliptischer Differentialgleichungen,

Teubner Studienbücher Mathematik, Stuttgart, 1986.

[12] Falgout R. D., Yang, U. M.: hypre: A library of High Performance Precondi-

tioners, lecture Notes in Computer Science, Springer, 2002.

38


[13] Kilian, S.; Münch, M.: A new generalized domain decomposition strategy

for the efficient parallel solution of the FDS-pressure equation, Part I:

Theory, Concept and Implementation. Technical Report ZR-09-19, Konrad-

Zuse-Zentrum für Informationstechnik Berlin, Juni 2009.

[14] Kilian, S.; Münch, M.: A new generalized domain decomposition strategy

for the efficient parallel solution of the FDS-pressure equation, Part !I:

Verification and Validation. Technical Report ZR-09-20, Konrad- Zuse-

Zentrum für Informationstechnik Berlin, Juni 2010.

[15] McDermott, R.: Notes on the Fire Dynamics Simulator (FDS) projection

scheme, Building and Fire Research laboratory, NIST Technical Note 1485,

2007.

[16] McGrattan, K. et. al.: Fire Dynamics Simulator (Version 5), Technical Refe-

rence Guide, Volumes 1-3, NIST Special Publication 1018-5.

[17] Saad, Y.: Iterative Methods for Sparse linear Systems, SIAM, Philadelphia,

2003.

WeItere InForMatIonen www.hhpberlin.de

hhpberlin ist das führende deutsche Ingenieurbüro für Brandschutz mit Sitz

in Berlin, München, Hamburg, Frankfurt und Braunschweig. Die im Jahr 2000

aus dem Büro Hosser, Hass und Partner hervorgegangene Firma entwickelt

weltweit Brandschutzkonzepte für nationale und internationale Bauprojekte.

Zu den Referenzen gehören u. a. die Münchner Allianz Arena, das Bundes-

kanzleramt, die Color line Arena in Hamburg, die Dalian Twin Towers, das

Pudong Museum in China sowie der neue Hauptstadt-Airport Berlin-Branden-

burg International. Die Kompetenz von hhpberlin reicht von der brandschutz-

gerechten Fachplanung über die Ausführung bis hin zur Qualitätssicherung

– sowohl im Neubau als auch bei der Bauerneuerung.

hhpberlin

Ingenieure für

Brandschutz GmbH

Hauptsitz:

Rotherstraße 19

10245 Berlin

Phone +49 (30) 895955 0

Fax +49 (30) 895955 9 101

www.hhpberlin.de

[email protected]

FDS: Numerische Optimierung und Kooperationen · Fassung 21. September 2009 autoren DR. RER. NAT....

Documents

Transcript of FDS: Numerische Optimierung und Kooperationen · Fassung 21. September 2009 autoren DR. RER. NAT....