Formatvorlage "UM-Vortragstitel"

Numerische Analyse eines 50000 rpm
High-Speed Turbogenerators
1
1
1
Hans-Christian Lahne , Volodymyr Bilyi , Tobias Drechsel ,
1
2
1
Jan Hiller , Oleg Moros , Dieter Gerling
1
Chair of Electrical Drives and Actuators, Universitaet der Bundeswehr Muenchen,
Werner-Heisenberg-Weg 39, 85577 Neubiberg, Germany
2
FEAAM GmbH, Werner-Heisenberg-Weg 39,
85579 Neubiberg
[email protected], [email protected]
Summary
This paper deals with the challenge of designing a high-speed induction machine with a rated speed of
50000 rpm for aircraft usage. Both the mechanical and the electrical aspects of the design are critical.
The simulation time is extensive due to the long period of time it takes to reach a steady-state
condition. In addition, the long simulation time is a result of the high-speed involved, together with its
small time steps. To understand these difficulties better, the following topics have been taken into
account:
•
•
•
•
•
The Comparison of Maxwell 2D and 3D Simulation Results: Time vs. Accuracy
Simulation-Specific Aspects Concerning Maxwell: Sampling in Maxwell
The Mechanical Analysis of the Model
Setting Up an Initial Current (among other things) to Reduce Simulation Time
The Usage of VBScript and Python
All results are based on analytical calculations and the FEM-based tools ANSYS Maxwell 2D, Maxwell
3D, and ANSYS Mechanical.
Keywords
High-Speed, High-Power, Aircraft Application, Turbo Generator, High Efficiency, Induction Machine,
Mechanical Stability, Simulation Time, VBScript, Python
ANSYS Conference &
33. CADFEM Users’ Meeting 2015
1
24. – 26. Juni 2015, Messe Bremen
1.
1.1
Vergleich der Simulationszeit und –genauigkeit von 2D- und 3D-Simulationen
Vergleich der 2D- und 3D-Simulationen
Da sich die Simulationszeiten von Maxwell 2D- und 3D-Berechnungen deutlich unterscheiden, soll
anhand der vorhandenen numerischen Ergebnisse eine Abschätzung der Abweichungen in Hinblick
auf die verwendete High-Speed Maschine vorgenommen werden. Ziel ist es, den Aufwand gegenüber
dem Nutzen in Relation zu stellen und damit zukünftige Simulationsaufgaben zeiteffizienter
durchführen zu können.
Bei der Modellierung in Maxwell 2D werden neben den Strangwiderständen und der Hauptinduktivität
die Streuinduktivitäten wie Nutstreuung und Oberwellenstreuung automatisch mit berücksichtigt, da es
2D-Effekte sind. Dagegen lassen sich die Stirnstreuung bzw. Wickelkopfstreuung nicht in 2D abbilden.
Der Wickelkopf fließt in die Simulation hierbei nur anhand seiner analytisch bestimmten
Wickelkopfinduktivtät mit ein [1]. Hier liegt der große zeitliche Unterschied zwischen Maxwell 2D und
3D. In 2D kommt nur der elektrische Einfluss der Wickelkopfinduktivität zum Tragen, wohingegen der
Einfluss der dadurch auftretenden Streufelder nur in einer numerisch-iterativen 3D-Simulation
Berücksichtigung findet.
Für die folgenden Betrachtungen wurde ein Turbogenerator ausgewählt, der im Rahmen eines
Luftfahrtprojekts mit einer Drehzahl von 50000 rpm bei einer maximalen Leistung von 700 kW
ausgelegt wird.
Analytische Betrachtungen hinsichtlich des Potentials verschiedener Legierungen und deren
Auswirkungen hinsichtlich der mechanischen Festigkeit finden sich in [2]. Demnach wurde der
Rotoraußendurchmesser rod unter Berücksichtigung eines Sicherheitsfaktors SF von 1,2 gewählt. Die
verwendeten Materialien und ggf. Blechdicke zeigt Tabelle 1.
Tabelle 1: In der Asynchronmaschine verwendete Materialien
Bereich
Statorwicklung
Rotorstäbe
Statoreisen
Rotoreisen
Material
Kupfer
BeCu
Vacoflux 48 - 0,1 mm
Vacodur S Plus - 0,35 mm
Beschreibungen der verwendeten Kobalt-Eisenlegierungen finden sich in [3] und [4]. Die Geometrie
der hierbei eingesetzten Asynchronmaschine mit Käfigläufer basiert auf Modellen aus [2] und [5],
wobei sowohl in letzterem, als auch in [6] relevante Designkriterien für High-Speed-Maschinen
beschrieben werden.
Im Folgenden werden drei Modelle erläutert, die in den Tabellen 2 bis 4 beschrieben werden. Modell 1
besitzt eine 36-2 Topologie, wohingegen Modell 2 eine 24-2 Konfiguration aufweist. Die erste Zahl
steht für die Anzahl der Statorzähne, wohingegen die Zweite die Anzahl der Pole angibt. Modell 3
besitzt ebenfalls eine 24-2 Konfiguration und ist eine Verbesserung von Modell 2, welches anhand von
Diskussionen gemäß [7] entstand. Auf die Geometrie wird aber an der Stelle nicht weiter
eingegangen. Für nähere Details wird auf [5] verwiesen.
Ein Problem der Asynchronmaschine ist der lange Einschwingvorgang bis zum Erreichen des
stationären Zustands [1]. Dies wirkt sich beträchtlich auf die Dauer der FEM-Simulationen aus.
Als Zeitschritt wurden 20 µs bei einer Endzeit von 0,12 s bis 0,15 s gewählt. Die 2D-Simulationszeit
des geschnittenen Modells betrug jeweils über 1 h, wohingegen die 3D-Modelle ca. 2 bis 4 Wochen
benötigten. Zur Sicherstellung gleicher Voraussetzungen wurden je Modellreihe gleiche Schlupfwerte
und Spannungen bei den Simulationen eingestellt. Die nachfolgenden Tabellen enthalten die
Ergebnisse der Simulationen.
Modell 1 (siehe Tabelle 2) stellt das erste Auslegungsergebnis dar, welches die hohe Drehzahl sowie
die Leistungsanforderungen erreichte. Das 3D-Modell weist ein knapp 2% geringeres Drehmoment als
die Referenz auf. Der Wirkungsgrad sank um knapp 1%. Die Ergebnisse des ursprünglichen RMXprtModells liegen - wie erwartet - näher an den Ergebnissen des 2D-Modells. Ein deutlich größerer
Unterschied ist bei den Drehmomentrippeln festzustellen.
ANSYS Conference &
33. CADFEM Users’ Meeting 2015
2
24. – 26. Juni 2015, Messe Bremen
Tabelle 2: Modell 1 mit 36-2 Konfiguration
Modell 1
Drehmoment [Nm]
Abweichung zum Maxwell2D-Wert [%]
Rippelmoment [Nm]
Wirkungsgrad
Simulationsdauer
RMXprt
134,9
-0,222
0,954
~1 s
Maxwell 2D
135,2
Referenz
7,8
0,958
2h
Maxwell 3D
132,9
-1,731
21,1
0,949
~ 30 d
Im Anschluss daran wurde ein verbessertes Modell mit einer 24-2 Konfiguration aufgebaut (siehe
Tabelle 3), welches insbesondere aus mechanischer Sicht Vorteile bietet.
Tabelle 3: Modell 2 mit 24-2 Konfiguration
Modell 2
Drehmoment [Nm]
Abweichung zum Maxwell2D-Wert [%]
Rippelmoment [Nm]
Wirkungsgrad
Simulationsdauer
RMXprt
141,7
5,646
0,959
~1 s
Maxwell 2D
133,7
Referenz
22,8
0,948
1 h 12 min
Maxwell 3D
135,4
1,256
21,3
0,946
~15 d
Bei Modell 3 (siehe Tabelle 4) unterscheiden sich die Drehmomentwerte von 2D zu 3D-Simulation nur
im geringen Prozentbereich. Insbesondere hier kann also die 3D Lösung sehr gut durch Maxwell 2D
abgebildet werden.
Tabelle 4: Modell 2 mit 24-2 Konfiguration
Modell 3
RMXprt
Maxwell 2D
Drehmoment [Nm]
Abweichung zum Maxwell2D-Wert [%]
Simulationsdauer
134,5
0,595
~1 s
133,7
Referenz
24 min
Maxwell 2D
(feines Mesh)
132,2
-1,135
1h 20min
Maxwell 3D
132,3
-1,058
21 d
Für die Ergebnisse aller betrachteten Modelle zeigt sich insgesamt nur eine minimale Abweichung des
Drehmoments. Dieser nahezu vernachlässigbare Unterschied rechtfertigt nicht den Mehraufwand an
Rechenzeit. Folglich kann als Ergebnis der durchgeführten Simulationen geschlossen werden, dass
bei guter Modellierung in Maxwell 2D die 3D-Simulationen im Designprozess weggelassen werden
können, insbesondere da keine Schrägung oder axiale Unsymmetrien in den Modellen auftauchen.
Nur beim späteren Bau eines Prototyps oder dem finalen Modell macht die detaillierte Berechnung
des 3D-Modells Sinn.
Ein Problem stellte hierbei die fehlende automatische Speicherfunktion von Maxwell 3D dar, welche
bei auftretenden Computerabstürzen oder Netzwerkabbrüchen zu einem Verlust der - bzw.
zutreffender formuliert - zu einem Verlust des Zugriffs auf die Simulationsdaten und einem Neubeginn
führte. Zwar konnten die bis zum Abbruch generierten Ergebnisse noch angeschaut werden, wodurch
klar wird, dass die Simulationsdaten an sich noch vorhanden waren, aber die neue Simulation begann
dennoch beim Anfangszeitpunkt. Als einzige Lösung kann hier in regelmäßigen Abständen ein „clean
stop“ durchgeführt und eine manuelle Speicherung vorgenommen werden. Falls dann ein Absturz
auftritt, fällt die Simulation nur auf diesen Speicherzeitpunkt zurück.
1.2
Signalabtastung in Maxwell
Um verwertbare Ergebnisse aus der FEM-Simulation zu bekommen, müssen ausreichend viele
Datenpunkte berücksichtigt werden. Dies ist vergleichbar mit der Abtastung von Signalen, welche zur
späteren exakten Rekonstruktion hinreichend oft abgetastet werden müssen.
Hierbei gilt das Nyquist-Abtasttheorem von Shannon, welches besagt [8], [9]:
ANSYS Conference &
33. CADFEM Users’ Meeting 2015
3
24. – 26. Juni 2015, Messe Bremen
1
 fs  2
Ts
fmax ;
fs
Abtastfrequenz; Ts
fmax
Höchste vorkommende Frequenz
Abtastintervall
(1)
Auf die zu berechnende Maschine angewendet heißt das wiederum, dass bei der Umdrehungszahl
1
von n  50000
eine Periode TPeriode
min
TPeriode

1
 1, 2
50000 1
60 s
103 s
(2)
dauert. Um hierbei eine ausreichende Abtastung zu gewährleisten, müssen also mehr als zwei
Datenpunkte pro Periode aufgenommen werden. Je höher die Abtastrate, desto zuverlässiger das
Ergebnis.
„Bei 5 Abtastungen pro Periode der höchsten Meßsignalfrequenz [...] beträgt der relative Abtastfehler
betragsmäßig noch 6,6%“ [12]. Falls dieser Abtastfehler Frel nur in einer Größenordnung von 1 bzw.
0,01% liegen soll, müssen dafür 12,8 bzw. 128 Abtastungen pro Periode stattfinden [12].
Insbesondere der Wert von 0,01% ist aber als rechnerischer, theoretischer Wert anzusehen.
Die folgende Abbildung 1 veranschaulicht diesen Zusammenhang zwischen relativem Abtastfehler
und dem Frequenzverhältnis
fmax
fs
. Ferner ist hier die Grenze des Shannonschen Abtasttheorems
1
eingetragen. Zur Berechnung des relativen Abtastfehlers eines Abtastkreises bei einer Frequenz f
dient folgende Formel [12]:
 f 
Frel  si    1
 fs 
(3)
Daraus ergibt sich eine Reihenentwicklung, die nach dem ersten Glied abgebrochen werden kann,
falls mehr als fünf Abtastungen pro Periode genutzt werden, da das Folgeglied weniger als 2% zum
Abtastfehler beiträgt. Hieraus kann dann die maximale Abtastfrequenz fmax als Funktion von fs
bestimmt werden [12]:
fmax
fs

1
6

Frel
(4)
Abbildung 1: Zusammenhang zwischen relativem Abtastfehler und dem Verhältnis von maximaler Frequenz zur
Abtastfrequenz [12]
1
Die Signalrekonstruktion aus den Abtastwerten erfolgt gemäß des Abtasttheorems mittels eines idealen Rechteckfilters.
Hierbei wird im Regelfall ein Abtast- und Haltekreis verwendet, bei dem der „abgetastete Wert bis zur nächsten Abtastung
beibehalten wird“ [12].
ANSYS Conference &
33. CADFEM Users’ Meeting 2015
4
24. – 26. Juni 2015, Messe Bremen
Gemäß dieser Überlegungen werden daher im Folgenden 128 Datenpunkte berücksichtigt. Das
notwendige Abtastintervall bzw. die Schrittweite Ts in Maxwell beträgt damit:
Ts 
1, 2
103 s
 9,375
128
106 s  1 105 s
Soll der Fehler maximal 0,1% betragen, ist ein Abtastintervall von rund 3
40 Abtastpunkten pro Periode entspricht.
(5)
105 s notwendig, was ca.
Im realen Betrieb weichen aber die Signalverläufe von den hier beschriebenen ab, da die auftretenden
Stromverläufe aus Grund- und Oberschwingungen bestehen, die jeweils unterschiedlich starke
Einflüsse auf das Maschinenverhalten haben. Daher müssen diese für jeden Fall genauer untersucht
werden. Dennoch dient die Betrachtung als Anhaltspunkt für jede Simulation.
Es ist anzumerken, dass Wellen eine räumliche Ausdehnung besitzen, wohingegen Schwingungen
zeitlich veränderliche Vorgänge beschreiben [1].
Die Ordnungszahlen des Signalspektrums entsprechen der folgenden Formel [1], wobei g 
6g  1  1, 5, 7, 11,13,...
:
(6)
1,000
MMF 24-2 und 36-2 [p. u.]
1,0
 24-2 Topologie
 36-2 Topologie
0,8
0,6
0,4
0,2
0,043
0,011 0,043
0,007
0,012
0,004
0,040
0,024
0,010 0,010
0,005
0,004
0,0
1 3 5 7 9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 45 47 49 51 53 55
Harmonic Order
Abbildung 2: MMF-Spektrum der Modelle 1 und 2 (normiert auf Grundwelle)
Die Abbildung 2 zeigt die MMF-Spektra der betrachteten Modelle 1 und 2. Hierbei ist festzuhalten,
dass sich beide Spektra relativ ähnlich sind. Das liegt daran, dass sich die Wicklungstopologien im
Wesentlichen nur durch die Statorzahnzahl unterscheiden und die Wicklungsfaktoren  der
Grundwellen damit sehr ähnlich sind ( 24 2  0,958; 362  0,956 ).
In Relation zur Amplitude der Grundwelle  ,1 lassen sich die Amplituden der Oberwellen (pro Einheit)
in Erweiterung einer Formel in [1] mit
 ,6g 1 
 ,6g 1  ,1
 ,1 6g  1
(7)
bestimmen, bei der die auftretenden Wicklungsfaktoren Berücksichtigung finden. Die kumulierten
Oberwellenanteile ergeben 18,2 %.
Im Vergleich dazu kann die „gesamte harmonische Verzerrung“ (THD, engl.: total harmonic distortion)
bestimmt werden, welche die Verzerrungen des Signals über die Zeit beschreibt [13], [14]. Diese kann
entweder in Relation zur Grundwelle oder zum Effektivwert des Signals gesetzt werden, wobei gemäß
Untersuchungen in [13] die erste Variante vorzuziehen ist. Es gilt damit [13]:
ANSYS Conference &
33. CADFEM Users’ Meeting 2015
5
24. – 26. Juni 2015, Messe Bremen

THD 
 Un
2
n 2
U1
 0, 064
(8)
Hierbei werden für U n die Spannungsamplituden der Oberwellen eingesetzt.
Aufgrund der gewählten Abtastrate sind die durchgeführten Simulationen demnach als hinreichend
genau einzustufen, da 64 bzw. zum Teil 128 Datenpunkte anstatt von nur mehr als zwei gemäß
Shannon berücksichtigt werden.
3.
Betrachtungen mechanischer Stabilität
Neben einer elektromagnetischen Simulation ist die Untersuchung der mechanischen Stabilität des
Rotors einer schnell rotierenden elektrischen Maschine von großer Bedeutung [10], [11].
Dazu kann eine gekoppelte numerische Simulation mittels Ansys Workbench durchgeführt werden.
Folgende Schritte werden dabei gemacht:








Erstellung eines 2D- bzw. 3D-Modells des Generators mit Ansys Maxwell
Import der Geometrie der untersuchten Maschine in Ansys Workbench
Definition der mechanischen Eigenschaften der Rotormaterialien
Definition der Kontaktflächen (Rotor – Stab – Kurzschlussring - Welle)
Erstellung einer Modellvernetzung
Definition der Rotation des Rotors
Durchführung der Simulation
Auswertung der Ergebnisse
Um die Simulationszeit zu reduzieren, ist es sinnvoll, nur ein Rotorsegment anstatt des gesamten
Rotors zu simulieren. Mit einer gekoppelten Simulation ist es möglich einen Auslegungsvorgang
zeitlich zu optimieren. Bei der Simulation geht es nur um eine statisch mechanische Analyse. Sowohl
die Wirkung der magnetischen Kräften als auch Temperatureffekte wurden dabei vernachlässigt.
3.1
Numerische Analyse mit Ansys Mechanical
In vielen Fällen lässt sich die Untersuchung der mechanischen Festigkeit von bestimmten Bauteilen
mit Hilfe von einem 2D FEM-Modell durchführen. In diversen wissenschaftlichen Veröffentlichungen
lassen sich Beispielmodelle finden. In diesen sind die Rotorstäbe nicht miteinander verbunden und die
Kurzschlussringe werden nicht betrachtet. Mit Hilfe solcher Annahmen kann die Simulationszeit
deutlich verringert werden. In der vorliegenden Arbeit wird für die genauere Untersuchung die
gesamte Rotorgeometrie von allen Komponenten mit berücksichtigt. Abbildung 3 zeigt das
vollständige 3D-Modell des Läufers von dem untersuchten Generator.
Abbildung 3: Komplettes 3D-Modell des Läufers auf der Welle
ANSYS Conference &
33. CADFEM Users’ Meeting 2015
6
24. – 26. Juni 2015, Messe Bremen
Der Rotor besteht aus einem Blechpaket in dem sich 18 Stäbe befinden, die mittels zweier
Kurzschlussringe verbunden sind. Das Ganze sitzt auf der Welle. Die mechanischen Eigenschaften
der für die Simulation verwendeten Materialien sind in Tabelle 5 zusammengefasst.
Tabelle 5: Materialeigenschaften für die statisch-mechanische Simulation
Name
Berylliumkupfer
(Rotorstäbe)
®
Vacodur S Plus
(Rotorblechpacket)
Wellenstahl CF53
(Welle)
Materialeigenschaften
Wert
3
Massendichte
8260 kg/m
E-Modul
Poissonzahl
115 GPa
0,3
Zug-Streckgrenze
689 MPa
Max. Zugfestigkeit
Dichte
800 MPa
3
8120 kg/m
E-Modul
Poissonzahl
250 GPa
0,3
Zug-Streckgrenze
750 MPa
Max. Zugfestigkeit
Dichte
1200 MPa
3
7800 kg/m
E-Modul
Poissonzahl
200 GPa
0,25
Zug-Streckgrenze
Max. Zugfestigkeit
430 MPa
830 MPa
Für die Reduzierung der Simulationsdauer wird nur 1/18 des Modells untersucht (siehe Abbildung 4).
Durch die Kennzeichnung der Ausschnittflächen mit dem Parameter „reibungsfreie Lagerung“ wird
eine Achsensymmetrie definiert.
a)
Abbildung 4: a) Modellausschnitt b) Generiertes Mesh
b)
Es wird angenommen, dass für die Fertigung der Rotorstäbe innerhalb des Rotorblechpakets ein
Druckgussverfahren verwendet wird. Die Kontakte zwischen Rotorstab und Rotorblechpaket und
zwischen Rotorblechpaket und Kurzschlussring sind als „Verbund“ und als „reibungsbehaftet“ mit
einem entsprechenden Reibungskoeffizienten von 0,2 definiert.
ANSYS Conference &
33. CADFEM Users’ Meeting 2015
7
24. – 26. Juni 2015, Messe Bremen
Die Simulation wird für eine Drehzahl von 50000 rpm durchgeführt. Abbildungen 5 bis 7 zeigen die
mechanische Spannung innerhalb des Rotors. Das Maximum der mechanischen Spannung im Rotor
beträgt ca. 708 MPa und liegt unter dem Grenzwert von 750 MPa. Die Abbildung 6 zeigt, dass die
Spannung innerhalb eines Rotorstabs maximal bei ca. 351 MPa liegt, sodass die Zug-Streckgrenze
von Berylliumkupfer nicht erreicht wird.
Abbildung 5: Mechanische Spannung innerhalb des Rotors
Abbildung 6: Mechanische Spannung im Rotorstab und Kurzschlussring
Abbildung 7 gibt Auskunft über die entstehende mechanische Spannung im Blechpaket des Rotors.
Der Wert von ca. 708 MPa wird sowohl durch die hohe Drehzahl als auch durch die Berücksichtigung
der Presspassung verursacht.
ANSYS Conference &
33. CADFEM Users’ Meeting 2015
8
24. – 26. Juni 2015, Messe Bremen
Abbildung 7: Mechanische Spannung im Rotorblechpaket
Die Gesamtverformung des Rotors wird in Abbildung 8 dargestellt. Die maximale Verformung von
0,09 mm tritt an den Kurzschlussringen auf. Die Welle erfährt eine geringere Deformation. Die
Gesamtverformung ist in dem Fall als unkritisch zu sehen.
Abbildung 8: Die Gesamtverformung des Rotors
Anhand der Simulationen wird deutlich, dass der maximale Wert der mechanischen Spannung für das
gewählte Material nicht überschritten wird. Die mechanische Stabilität des Rotors ist damit
gewährleistet.
3.2
Festigkeitsnachweis der Welle mittels Analyse einer Presspassung
Neben den Critical Speeds gilt es im Rahmen der mechanischen Festigkeitsanalyse u.a. auch die
Festigkeit der Presspassung zu untersuchen. Um den Rotor auf die Nenndrehzahl von 50000 rpm zu
beschleunigen, muss die Welle das jeweilige Drehmoment übertragen können, ohne dass ein
Durchrutschen stattfindet. Der Kontaktdruck zwischen der Welle und dem Blechpaket muss stets so
groß sein, dass der Zustand der Kontaktfläche als „haftend“ bezeichnet werden kann.
In der folgenden Simulation in der Ansys Workbench wird der Rotor als Segment-Modell dargestellt.
Das wiederum hat den Vorteil einer geringeren Berechnungszeit. Die Welle hat gegenüber dem
Blechpacket ein Übermaß von 50 μm. Bei einem angenommenen Reibungskoeffizienten von 0,2
ANSYS Conference &
33. CADFEM Users’ Meeting 2015
9
24. – 26. Juni 2015, Messe Bremen
zwischen Wellenstahl und dem aus Vacodur S Plus gefertigten Blechpacket wird die Wirkung des
Kontaktdrucks bei Nenndrehzahl untersucht.
Dazu wird vereinfachend angenommen, dass die Welle an der einen Seite im Mittelpunkt fest gelagert
wird. An der anderen Seite ist die Welle so gelagert, dass eine Verformung in z-Richtung möglich ist.
Eine Ausdehnung in z-Richtung erfolgt in der Praxis beispielsweise durch thermische Verformung. Bei
der Spannungsanalyse ist darauf zu achten die Vernetzung so zu verfeinern, dass man pro
Spannungsfall mindestens drei Knoten erhält.
Eine Simulation der Kontaktdruckverteilung liefert die Ergebnisse auf Abbildung 9:
Abbildung 9: Kontaktdruckverteilung in Ansys Mechanical
Betrachtet wird nur die Kontaktfläche zwischen der Welle und dem Blechpaket. Abbildung 9 zeigt
dabei eine Draufsicht der Welle. Dabei fällt auf, dass der Kontaktdruck in der Mitte der Welle bei ca.
99 MPa liegt und erst an den Seiten des Blechpakets abnimmt. An den beiden Randstellen beträgt der
Kontaktdruck ca. 2 MPa.
Der Zustand der Kontaktfläche wird bei Ansys Mechanical mit einem Kontakt-Tool dargestellt und
kann aufgrund des Kontaktdrucks für die gesamte Fläche als „nicht haftend“ bezeichnet werden. Der
Kontaktdruck reicht nicht aus, um ein Durchrutschen der Welle beim Betrieb der Maschine zu
verhindern und um die mechanischen Drehmomente zu übertragen. Hieraus resultiert, dass ein
anderes Verfahren benutzt werden muss.
4.
Möglichkeiten zur Reduktion der Simulationszeit bzw. zum schnelleren Erreichen eines
stationären Zustands
Im Folgenden sollen mögliche Ansätze zum schnelleren Erreichen eines eingeschwungenen Zustands
und Lösungen zur Minimierung der Simulationszeit vorgestellt werden. Die Vorgehensweise und
Ergebnisse sollen in diesem Kapitel erläutert werden. Im Einzelnen sind das die Definition von
Anfangsstromwerten, die Vorgabe eines Rotorpositionswinkels, eine Cosinus- und eine ExponentialFunktion und schließlich eine Optimierung der Zeitschrittweite.
4.1 Definition von Anfangsstromwerten
Eine denkbare Möglichkeit zur Reduktion der Simulationszeit ist die Vorgabe von Anfangsstromwerten
„Initial current“, falls eine Spannungsspeisung eingestellt wird. Diese können sowohl für eine 3D-, als
auch für eine 2D-Simulation aus einer vorhergehenden 2D-Simulation übernommen werden, bei der
der eingeschwungene Zustand bereits erreicht wurde. Das Verfahren zum Auslesen der
Anfangsstromwerte wird in Abbildung 10 dargestellt. Der gewählte Zeitpunkt wird durch die grüne
Linie markiert. Hierbei werden die drei Phasenströme zu einem Nulldurchgang der Spannung von
Phase A aufgenommen und in die neue Simulation übernommen.
ANSYS Conference &
33. CADFEM Users’ Meeting 2015
10
24. – 26. Juni 2015, Messe Bremen
Abbildung 10: Auslesen der Anfangsstromwerte aus Maxwell
Im direkten Vergleich konnte festgestellt werden, dass sich ein größerer Überschwinger im
Drehmomentverlauf zu Beginn einstellte (siehe Abbildung 11), der um den Faktor 20 größer war als in
der Referenzsimulation (siehe Abbildung 13). Zudem reduzierten sich weder die Zeit bis zum
Erreichen des stationären Zustands, noch die Gesamtsimulationszeit. Als weiterer negativer Punkt
spalteten sich die Stromverläufe – wie in Abbildung 17 zu sehen - weiter auf.
Abbildung 11: Drehmomentverlauf bei Eingabe von Anfangsstromwerten
Dementsprechend ist dieser Ansatz bisher nicht zielführend und bedarf weiterer Überprüfung, da die
Theorie eine Verbesserung in Aussicht stellt.
4.2 Vorgabe eines Rotorstartwinkels
Im Vergleich dazu wurde ein anderer Ansatz gewählt, bei dem als Erweiterung des Ansatzes in
Kapitel 4.1 der Zeitpunkt des stationären Zustands ausgelesen und mit als Startwert für die
Spannungsspeisung in allen drei Phasen übernommen wurde. Unter Zuhilfenahme dieses Zeitpunkts
wurde ein neuer Rotorstartwinkel berechnet. Zur Bestimmung des Rotorstartwinkels wurde die
Drehzahl mit dem Zeitwert des stationären Zustands multipliziert, der ganze Zahlenanteil entfernt und
die Dezimalstellen mit 360° multipliziert. Bei diesem Vorgehen war keine Veränderung zum vorherigen
ANSYS Conference &
33. CADFEM Users’ Meeting 2015
11
24. – 26. Juni 2015, Messe Bremen
Kapitel feststellbar. Folglich sind für Kapitel 4.1 und Kapitel 4.2 nur gemeinsame Abbildungen
aufgeführt.
4.3 Verwendung einer Cos-Funktion
Als nächste Maßnahme wurde anstatt der üblicherweise verwendeten Sin-Speisung eine
Spannungsspeisung mit Cosinus vorgenommen. Dieser Versuch beruht auf der Überlegung, dass die
Cosinus-Funktion beim Maximum startet und kleiner wird, wohingegen die Sinus-Funktion vom
Minimum auf das Maximum steigt, wodurch gerade bei anfangs stark aufspreizenden Strömen eine
Verbesserung erzielt werden könnte. Die Auswirkungen bewegten sich jedoch in vernachlässigbaren
Rahmen. Daher werden hier keine Beispielbilder präsentiert.
4.4 Eingabe einer Exponentialfunktion
Als weitere Idee zur Optimierung wurde bei der Spannungsspeisung eine Exponentialfunktion
 time
3
( 1  e 3 10 [7]) hinzugefügt, die den Einschwingvorgang beschleunigen sollte. Ergebnisse zeigen
signifikant reduzierte Drehmomentrippel während des Einschwingens (siehe Abbildung 19), die Zeit
bis zum Erreichen des stationären Zustands sinkt von ca. 70 ms auf ca. 60 ms, wobei die Ströme
verglichen mit der Referenz einen deutlich reduzierten Ausschlag aufweisen. Gibt es bei der Referenz
noch Amplituden bis fast 12500 A, weisen die Ströme der optimierten Simulation Maximalwerte von
unter 2500 A auf (vergleiche Abbildung 19 und Abbildung 20). Es müssen aber durch diese Funktion
Simulationsungenauigkeiten in Kauf genommen werden. Ferner scheint ein nicht unerheblicher Rippel
bei den Verlusten (Abbildung 21) zu entstehen, deren Mittelwert sich jedoch wiederum kaum vom
Ursprungswert unterscheidet. Ein Vergleich der Verlustamplituden ergibt allerdings, dass die Amplituden von zu Beginn 300 kW (Abbildung 15), auf knapp unter 16 kW sinken (Abbildung 21), was ein
bedeutender Vorteil ist. Ferner lässt sich erneut klar erkennen, dass der eingeschwungene Zustand
früher erreicht wird. Daher wird eine positive Bewertung gegeben. Insgesamt stellt dieser Ansatz eine
sinnvolle Möglichkeit dar, um bei zeitintensiven Simulationen Optimierungen vorzunehmen.
4.5 Optimierung der Zeitschrittweite
Die letzte Optimierungsmaßnahme zielt hauptsächlich auf die Reduktion der gesamten
Simulationszeit ab. Hierzu wurden variable Zeitschrittweiten, anstatt der sonst üblichen konstanten
Zeitschritte genutzt. Im Regelfall wird ein fixer Zeitschritt gemäß der Überlegungen aus Kapitel 1.2
gewählt. Stattdessen wird eine Funktion eingegeben, die am Anfang größere Schritte vorsieht und
dann gegen einen konstanten Wert konvergiert (siehe Abbildung 12). Allgemein formuliert lautet diese
empirisch ermittelte Funktion:
f (time )  9
time _ step
0,5


 time  
1  cos 
  time _ step 
 stop _ time 


(9)
Die Laufvariable ist die Zeit „time“, die Zeitschritte entsprechen „time_step“ und das Simulationsende
ist durch „stop_time“ definiert. Zusätzlich wurde eine Funktion
time _step  if  time  stop _time2, formel," kons tante" 
(10)
implementiert, die bei einer Zeit kleiner der „stop_time2“, die oben genannte Formel anwendet und für
einen erwarteten Zeitbereich des eingeschwungenen Zustands auf den unter „konstante“ definierten
Wert wechselt.
Als Ergebnis weisen die Stromverläufe in Abbildung 23 größere Amplituden als die Referenz in
Abbildung 14 auf. Ferner verschlechterten sich der Drehmomentverlauf (siehe Abbildung 22) und die
reale Zeit bis zum Erreichen des stationären Zustands geringfügig. Die zuletzt genannte Zeit stieg von
ca. 70 ms auf ca. 90 ms. Dagegen konnte aber die gesamte Simulationszeit deutlich reduziert werden.
Dieser Umstand ist hier der entscheidende, innovative Aspekt dieser Maßnahme. Die
Simulationsdauer sank auf 1/3 des ursprünglichen Werts. Gerade bei den in den Tabellen 2 bis 4
gezeigten extremen Zeiten der 3D-Simulationen muss sich diese Optimierungsmaßnahme sehr positiv
auswirken.
ANSYS Conference &
33. CADFEM Users’ Meeting 2015
12
24. – 26. Juni 2015, Messe Bremen
Zeitschrittweite


Abbildung 12: Optimierungsfunktion für die Zeitschritte
Zeitschritt [s]
4.6 Auswahl eines Ansatzes und Gesamtbewertung
Werden nun die Vorgehensweisen von Kapitel 4.4 und 4.5 kombiniert, lassen sich Simulationen
deutlich effizienter durchführen und merklich verkürzen. Die deutlichen Unterschiede verdeutlicht der
direkte Vergleich in Abbildung 25. Sowohl die Zeit bis zum Erreichen des stationären Zustands, als
auch die Gesamtsimulationsdauer sinken. Unter anderem können anhand der Auswertung von
betreffenden Verläufen vor Abschluss der eigentlichen Simulation Abschätzungen getroffen werden, in
welchem Wertebereich die stationäre Lösung liegen wird. Insbesondere für 3D-Simulationen stellt
diese Vorgehensweise einen entscheidenden Fortschritt dar.
Die Bewertungsmatrix (Tabelle 6) fasst die Ergebnisse aller Methoden zusammen und gibt eine
Aussage über die Eignung für zukünftige Simulationen. Hierbei ist die Gesamtdauer der Simulation
das ausschlaggebende Hauptkriterium. Eine positive Bewertung dieser Eigenschaft zieht automatisch
eine positive Eignung nach sich. Eine „0“ im Feld „Gesamtdauer der Simulation“ überträgt die
Bewertung der Eignung auf das Feld „Einschwingzeit“, welches die Zeit bis zum Erreichen eines
stationären Zustands beschreibt. Die Einzelkriterien sind indifferent hinsichtlich der Eignung und
dienen ausschließlich der Übersicht.
Im Anschluss an die Bewertungsmatrix folgen repräsentative Simulationsergebnisse.
Tabelle 6: Bewertungsmatrix hinsichtlich der vorgestellten Methoden zur Optimierung von Simulationen
Einzel-Kriterien
Eigenschaften
DrehmomentStromverlauf
verlauf
Optimierungsziele
Verlustkennlinien
Kapitel
Ansatz
4.1
Anfangsstromwert
-
-
0
4.2
Vorgabe
Rotorpositionswinkel
0
0
4.3
Cosinus-Funktion
0
4.4
Exponentialfunktion
4.5
4.6
Einschwingzeit
Gesamtdauer
der Simulation
0
0
-
0
0
0
0
+
+
+
+
0
Variabler Zeitschritt
0
-
0
0
+
Variabler Zeitschritt +
Exponentialfunktion
+
+
+
0
+
0
Legende:
+ = Verbesserung
- = Verschlechterung
ANSYS Conference &
33. CADFEM Users’ Meeting 2015
13
24. – 26. Juni 2015, Messe Bremen
0 = keine oder wenig Veränderung
Eignung
0
+
+
+
Referenzsimulation:
Abbildung 13: Drehmomentverlauf der Referenzsimulation
Abbildung 14: Stromverlauf der Referenzsimulation
Abbildung 15: Verlustverläufe der Referenzsimulation
ANSYS Conference &
33. CADFEM Users’ Meeting 2015
14
24. – 26. Juni 2015, Messe Bremen
Mit Anfangsstromwert und Rotorpositionswinkel:
Abbildung 16: Drehmomentverlauf mit Anfangsstromwert und Rotorpositionswinkel
Abbildung 17: Stromverlauf mit Anfangsstromwert und Rotorpositionswinkel
Abbildung 18: Verlustverläufe mit Anfangsstromwert und Rotorpositionswinkel
ANSYS Conference &
33. CADFEM Users’ Meeting 2015
15
24. – 26. Juni 2015, Messe Bremen
Mit Exponentialfunktion:
Abbildung 19: Drehmomentverlauf mit Exponentialfunktion
Abbildung 20: Stromverlauf mit Exponentialfunktion
Abbildung 21: Verlustverläufe mit Exponentialfunktion
ANSYS Conference &
33. CADFEM Users’ Meeting 2015
16
24. – 26. Juni 2015, Messe Bremen
Mit empirisch-ermittelter Funktion:
Abbildung 22: Drehmoment mit empirisch-ermittelter Funktion
Abbildung 23: Stromverlauf mit empirisch-ermittelter Funktion
Abbildung 24: Verlustverläufe mit empirisch-ermittelter Funktion
ANSYS Conference &
33. CADFEM Users’ Meeting 2015
17
24. – 26. Juni 2015, Messe Bremen
Direkter Vergleich der Simulationsergebnisse:
Abbildung 25: Direkter Vergleich von Referenz (links) zur Kombination der Exponential- und empirischen Funktion
(rechts)
5.
Visual Basic, VBScript und Python
In jeder Auslegungsphase müssen in den Programmen Maxwell 2D und 3D jeweils viele geometrische
Parameter wie Statoraußen- und Innendurchmesser, Jochbreite, Luftspalt und elektrische Parameter
wie Strangspannung, Speisefrequenz, u.v.m. eingestellt und optimiert werden. Um die Eingabe zu
erleichtern und viel Zeit zu sparen, wird daher im Folgenden für beide Programme eine
Parametrierung mit Hilfe von Visual Basic Script [15], [16] vorgestellt.
„Visual Basic Scripting Edition, kurz VBScript (VBS), ist eine Skriptsprache, welche für die WebProgrammierung sowohl auf Client- als auch auf Serverseite eingesetzt werden kann“ [17]. Die Syntax
von VBScript orientiert sich an Visual Basic und Visual Basic for Applications (VBA) [18].
Beide in diesem Bericht vorgestellten Skripte basieren auf einem früher erstellten Visual Basic Script,
welches zur Erstellung kompletter Modelle in Maxwell diente, deren Erstellung über RMxprt nicht
möglich ist [19]. Der Hauptunterschied ist, dass die beiden neuen Skripte zum einen bestehende
Modelle parametrisieren, anstatt jedes Mal ein neues Modell zu generieren und zum anderen
bestehende Maschinenparameter auslesen und bei der darauffolgenden Parametrisierung nutzen
können. Eine Neuerstellung eines gesamten Grundwellenmodells mittels VBS ist in Maxwell nicht
zielführend, da einfache Grundwellenmaschinen bereits gut durch RMXprt erstellbar sind. VBS kann
wiederum genutzt werden, um die Modelle von Oberwellenmaschinen in Maxwell aufzubauen.
Alternativ kann aber auch die Geometrie dieser Art von Maschinen in RMXprt erstellt und im Zuge der
Parametrisierung mit einer passenden Wicklungstopologie versehen werden. Im Vergleich von beiden
Ansätzen stellt diese Vorgehensweise den kürzeren Weg dar.
Abbildung 26 zeigt einen Teil der realisierten Parametrisierung eines bestehenden Maxwell 2DModells mit Hilfe von VBS. Damit ist eine einfache Änderung von Geometrie- oder Simulationsdaten
direkt in Maxwell möglich.
ANSYS Conference &
33. CADFEM Users’ Meeting 2015
18
24. – 26. Juni 2015, Messe Bremen
Abbildung 26: Parametrisierungsbeispiel mit VBS
Im Folgenden soll nun Python erläutert werden. Python ist eine Programmiersprache, die von Beginn
an objektorientiertes Programmieren ermöglicht, aber gleichermaßen für prozedurale Programmierung
bestens geeignet ist. Der Baustein Tkinter umfasst ein Werkzeug, mit dem grafische
Benutzeroberflächen plattformübergreifend generiert werden können. Das Programm Maxwell von
Ansys bietet zunehmend mehr Unterstützung für die Skriptsprache IronPython an, welches eine
Unterform von Python darstellt [20]. So können seit einigen Maxwell-Versionen Prozesse als Skript
aufgenommen werden.
Um das volle Potential der Skripte ausnutzen zu können, müssen sie für den jeweiligen Einsatz durch
Nutzereingaben konfigurierbar sein. Bisher wurde dies meist durch einzelne, aufeinanderfolgende
Eingabeboxen realisiert, welche zur Laufzeit des Skripts aufgerufen wurden.
Abbildung 27: Oberfläche des Python-Programms
Dabei war die schwierige Handhabung, fehlende Überprüfung der Eingaben und eine fehlende
Möglichkeit zum Abspeichern der Ergebnisse sehr problematisch. Die bisherige Verwendung von VBSkript ist in Zukunft nicht mehr möglich, da Microsoft die Unterstützung eingestellt hat.
Abbildung 27 zeigt einen Ausschnitt der neuen Benutzeroberfläche. Die Eingaben durch den Nutzer
werden bereits vor dem eigentlichen Erstellen des Python Skripts auf Fehler überprüft. Sollte ein
Fehler auftreten, erhält der Nutzer eine individualisierte Fehlermeldung. Weitere Lösungsvorschläge
und ausführliche Beschreibungen befinden sich in der Hilfe des Programms.
Nach dem Ausfüllen des Formulars lässt sich ein Python Skript erstellen und abspeichern. Dieses
Skript kann nun komplett ohne zusätzliche Nutzereingaben von Maxwell ausgeführt und zur
Parametrisierung eines Modells verwendet werden. Alternativ kann das Programm auch direkt
ANSYS Conference &
33. CADFEM Users’ Meeting 2015
19
24. – 26. Juni 2015, Messe Bremen
Maxwell starten und das Skript ausführen, dadurch ergibt sich die Möglichkeit verschiedene Versionen
des Skripts schnell und effektiv zu testen.
Falls Änderungen an den Einstellungen erforderlich werden, lassen sich alle Informationen eines
bereits erstellten Skripts durch das Programm einlesen und nachträglich verändern.
Das Programm leitet unerfahrene Nutzer durch die Erstellung einer skriptbasierten
Modellparametrisierung, gleichzeitig lässt es alle Möglichkeiten offen, das Skript eigenständig zu
erweitern. Die gesamte Oberfläche und Programmlogik sind in Java umgesetzt und dadurch
plattformunabhängig.
6.
[1]
[2]
[3]
[4]
[5]
[6]
[7]
[8]
[9]
[10]
[11]
[12]
[13]
[14]
[15]
[16]
[17]
[18]
[19]
[20]
Quellen
Gerling, D.: “Electrical Machines - Mathematical Fundamentals of Machine Topologies”,
Springer Verlag Berlin Heidelberg, Germany, 2014.
Lahne, H.C., Gerling, D.: “Investigation of High-Performance Materials in Design of a 50000
rpm High-Speed Induction Generator for use in Aircraft Applications”, Workshop on Aircraft
System Technologies (AST-2015), Hamburg, Germany, 24.-25.02.2015.
Volbers, N., Gerster, J.: High Saturation, High Strength Iron-Cobalt Alloy for Electrical
Machines, Proceedings of the INDUCTICA, CWIEME Berlin 2012.
Fohr, F., Volbers, N.: „Hochsättigende Eisen-Kobalt-Legierungen mit erhöhter Festigkeit“,
FEMAG Anwendertreffen, 09.11.2012.
Lahne, H.C., Gerling, D.: “Comparison of State-of-the-Art High-Speed High-Power
Machines”, 41st Annual Conference of the IEEE Industrial Electronics Society IECON2015,
Yokohama, Japan, 09.-12.11.2015 (unpublished).
Lahne, H.C., Gerling, D.: “Considerations in designing a 50000 rpm high-speed high-power
machine”, EPE’15-ECCE Europe, Geneva, Switzerland, 08.-10.09.2015 (in press).
Rudnyi, E., Bachinski-Pinhal, D., Schönborn, K., Bock, U.: CADFEM Support,
anwenderbezogene Kommunikationen, CADFEM GmbH, 2013-2015.
Grote, K. H., Feldhusen, J.: „Dubbel – Taschenbuch für den Maschinenbau“, 23. Auflage,
Springer-Verlag Berlin Heidelberg, 2011.
Höher, P. A.: „Grundlagen der Informationsübertragung“, Von der Theorie zur
Mobilfunkanwendung, zweite überarbeitete und erweiterte Auflage, Springer Vieweg Verlag,
Wiesbaden, 2013.
Binder, A., Schneider, T., Klohr, M.: ”Fixation of buried and surface-mounted magnets in
high-speed permanent-magnet synchronous machines”, IEEE Transactions on Industry
Applications, vol. 42, No. 4, July/August 2006, pp. 1031-1037.
Lovelace, E. C., Jahns, T. M., Keim, T. A., Lang, J. H.: “Mechanical design considerations
for conventionally laminated, high-speed, interior PM synchronous machine rotors”, IEEE
Transactions on Industry Applications, vol. 40, No. 3, May/June 2004, pp. 806-812.
Tränkler, H.-R.: „Taschenbuch der Messtechnik“, 5. Auflage, Oldenbourg Verlag, München,
2005.
Shmilovitz, D.: “On the definition of total harmonic distortion and its effect on measurement
interpretation”, IEEE Transactions on Power Delivery, vol. 20, no. 1, pp. 526-528, 14.02.2005.
“IEEE Standard Definitions for the Measurement of Electric Power Quantities Under
Sinusoidal, Nonsinusoidal, Balanced, or Unbalanced Conditions – Redline”, 12.07.2011.
Aitken, P.: „JavaScript und VBScript,“ Thomson Verlag, Bonn 1997.
Schwarz, R., Malluf, I.: „VBScript“, München, 1997.
Avci, O., Trittmann, R., Mellis, W.: „Web-Programmierung, Softwareentwicklung mit
Internet-Technologien - Grundlagen, Auswahl, Einsatz“ - XHTML & HTML, CSS, XML,
JavaScript, VBScript, PHP, ASP, Java, Vieweg Verlag, Wiesbaden, 2003.
Clark, S., De Donatis, A., Kingsley-Hughges, A., Kingsley-Hughes, K., Matsik, B.,
Nelson, E., Prussak, P., Read, D., Thomsen, C., Updegrave, S., Wilton, P.: “VBScript –
Programmer’ s reference”, Birmingham 1999.
Moros, O.: „Visual Basic Script zur Erstellung von Maxwell Maschinen in 24-10 bzw. 36-10
Konfiguration“, private Kommunikation, FEAAM GmbH, UniBw München, 2012.
http://www.pythonmania.de/article/warum.html
ANSYS Conference &
33. CADFEM Users’ Meeting 2015
20
24. – 26. Juni 2015, Messe Bremen