Computersimulation der induktiven Erwärmung

Die induktive Erwärmung kann sehr gut mithilfe der Computertechnik simuliert werden. Die Simulationen verbessern physikalische Prozesse der induktiven Erwärmung zu verstehen und sie gewähren wertvolle Informationen über Größen, die experimentell oft nur sehr schwierig messbar sind. Wir können so in das Werkstück hineinschauen und feststellen, wo die Wärme entsteht und wohin sie geleitet wird. Im Unterschied zum Experiment lassen sich mithilfe der Simulation ganze Abhängigkeiten von Eingangsgrößen feststellen, was bei dem Entwurf und der Optimierung genutzt wird. Die häufigsten Aufgaben sind: Festsetzung der Zeit, die zur Durchwärmung des Werkstück gebraucht wird, Ausrechnung des Temperaturprofils Kern-Oberfläche, Optimierung der Spulenform, Festsetzung der Erwärmungseffizienz. Die Simulationen helfen uns nicht nur bei dem Entwurf einer neuen Erwärmungsanlage, sondern auch bei der Feststellung der optimalen Einstellung des Erwärmers im Betrieb. Wir nutzen einige Berechnungsprogramme je nach der Komplexität der Aufgabe. Für die Simulation unserer Erwärmer entwickeln wir ein eigenes Berechnungsprogramm.

Schicken Sie uns Ihre Anfrage über die Simulation der Erwärmung an indukce@roboterm.cz oder per unser Anfrageformular:

Beispiele der Simulation der induktiven Erwärmung

Simulation der induktiven Erwärmung eines vierkantigen Blockes
Simulation der induktiven Erwärmung eines Brennteils aus Stahl
Simulation der induktiven Erwärmung des Schrankes für Aufschrumpfen der Lager

Wie läuft die Simulation ab

Aus Sicht des Simulationsingenieurs kann der Prozess in drei Phasen geteilt werden: Datenvorbereitung (preprocessing), eigene Berechnung und Verarbeitung der Ergebnisse (postprocessing). Die Datenvorbereitung besteht in der Schaffung eines geometrischen Modells. Dieses wird durch ein Netz in feine Elemente (gewöhnlich Dreiecke in 2D Aufgabe, Tetraeder in 3D Aufgabe) geteilt. Den Zonen müssen Materialeigenschaften (z. B. Wärmeleitfähigkeit), Anfangsbedingungen (z.B. Anfangstemperatur) und Randbedingungen (z.B. Strom in der Spule, Temperaturverlust durch Wärmestrahlung) zugeordnet werden. Der Simulationsingenieur wählt das passende Berechnungsmodell und stellt seine Parameter ein (z.B. Erwärmungszeit, Zeitschritt). Die Simulation der induktiven Erwärmung ist ein kombiniertes Problem, in dem das elektromagnetische Feld mithilfe der harmonischen Analyse und das Wärmefeld als Übergangsprozess gelöst werden. Die Entstehung der Wärme im Material entspricht dabei den Jouleschen Verlusten. Das Berechnungsprogramm beschreibt jedes Element mit Differenzialgleichungen (Maxwell-Gleichungen, Wärmeleitungsgleichungen) und löst dieses Gleichungssystem nummerisch in jedem Zeitschritt. Der Output ist ein zeitabhängiges elektromagnetisches und Wärmefeld. Das Ergebnis wird schließlich vom Berechnungsingenieur im Postprocessor in eine geeignete grafische Form aufbereitet (Graph, Animation).

Genauigkeit der Simulation

Oft werden wir gefragt, in wieweit das Ergebnis der Simulation der Wirklichkeit entspricht. Die Genauigkeit des Ergebnis beeinflussen mehrere Faktoren: Wahl des Berechnungsmodells, Vereinfachung der Geometrie, Netzdichte, Größe des Zeitschritts, Kenntnis der Materialeigenschaften und der Rand- und Anfangsbedingungen. Der Simulationsingenieur muss das geeignete Berechnungsmodell wählen. Wenn der magnetische Stahl auf die Schmiedetemperatur erwärmt wird, muss die Software mit einem nichtlinearem Algorithmus benutzt werden, der mit dem starken Rückgang der Permeabilität gegen die Curie-Temperatur rechnet. Ein Nachteil von solch einem Algorithmus ist jedoch, dass er erheblich die Laufzeit der Berechnung verlängert. Im Falle der hohen Anzahl der Elemente (3D Aufgaben, schwierigere 2D Aufgaben) stoßen wir bald an die Grenzen der üblichen Computertechnik. Um die Berechnung zu beschleunigen ist nötig die Aufgabe zu vereinfachen. 3D Aufgaben lassen sich oft in 2D Aufgaben umwandeln. In einigen Fällen, wie. z.B. bei der Optimierung der Spulenform, ist nicht nötig, den ganzen Übergangsprozess zu berechnen, sondern hier reicht nur die harmonische Analyse (die Temperaturen werden nicht berechnet, sondern nur Joulesche Wärme). Der Simulationsingenieur muss entscheiden, welche Vereinfachung er akzeptieren kann, ohne dass das Ergebnis deutlich beeinflusst wird. Das Netz muss sehr genau entwickelt werden. Das Netz muss ausreichend fein sein, jedoch nicht zu viel, damit die Laufzeit nicht zu lang wäre. Deshalb wird das Netz mit verschiedener Dichte entworfen, wobei das dichteste Netz in den Stellen mit den größten Veränderungen der berechneten Größen sein muss. Besonders fein muss das Netz in der Oberflächenschicht des Werkstückes im Fall der kleinen Eindringtiefe wie z.B. beim Stahl im magnetischen Stand. Wenn der Simulationsingenieur ein geeignetes Simulationsmodell schafft, hängt die Genauigkeit des Ergebnisses nur von der Genauigkeit der Eintrittsdaten (Materialeigenschaften, Rand- und Anfangsbedingungen) ab. Gerade die Nichtkenntnis der genauen Eintrittsdaten ist die Hauptursache der oft diskutierten Differenz zwischen der Simulation und der Wirklichkeit. Weil die Feststellung der Eintrittsdaten (besonders Materialeigenschaften und Emissionsgrad) im ganzen Temperaturbereich aufwendig ist, wird für die Berechnung eine Schätzung benutzt. Die genauen absoluten Ergebnisse können erst nach der Einstellung der Eintrittsdaten nach dem Experiment erreicht werden. Oft sind für die richtigen Ergebnisse nicht die absoluten Werte wesentlich, sondern vor allem die Abhängigkeiten. Die können mithilfe der Simulationen sehr zuverlässig und auch ohne Kenntnis der genauen Eintrittsdaten ausgewertet werden.

Abhängigkeit der Werkstoffeigenschaften von der Temperatur

Die Werkstoffeigenschaften können stark von der Temperatur abhängig sein. Wenn das Ergebnis der Simulation genauer der Wirklichkeit entsprechen soll, muss die Simulation mit diesem Fakt rechnen. Bei Stählen sind die Temperaturcharakteristiken stark von ihren chemischen Zusammensetzungen abhängig. Die Graphen unten zeigen, wie die Materialeigenschaften des Konstruktionsstahls mit dem Kohlenstoffgehalt von 0,2 % von Temperatur abhängen können. Aus den Diagrammen ergibt sich, dass um 760°C (Curie-Temperatur) zur Phasenänderung kommt. Bei dieser Temperatur verliert das Material ihre ferromagnetischen Eigenschaften und die relative Permeabilität sinkt heftig auf 1. Die Permeabilität der ferromagnetischen Stoffe lässt sich bei Berechnungen nur schwer beschreiben, weil sie außer der Temperatur auch von der Intensität des magnetischen Feldes abhängig ist. Der Graf der Permeabilität zeigt hier nur ihren Rückgang, die absoluten Werte würde man von der Hysteresekurve des Materials erfahren. Die nächste problematische Eigenschaft ist der Emissionsgrad, die die Menge der Wärmeverluste angibt, die durch Oberfläche ausgestrahlt wird. Der Emissionsgrad ist von der Materialoberfläche abhängig. Die Oberfläche des Stahls oxidiert bei hohen Temperaturen, was den Emissionsgrad erhöht.


Bild 1 - Spezifische Wärmekapazität	Bild 2 – Wärmeleitfähigkeit

Bild 3 - Resistivität	Bild 4 – Relative Permeabilität

Theorie der 1D Erwärmung

Um die Aufgabe zu vereinfachen, beschreiben wir ein drehsymmetrisches 1D Problem, das der statischen Induktionserwärmung eines langen zylindrischen Werkstückes in der langen Spule entspricht. Setzen wir voraus, dass die Materialeigenschaften konstant bleiben. Wir lösen also den Verlauf der physikalischen Größen nur in der Radialrichtung. Gehen wir von der Ampereschen (1) und Faradayschen (2) Gesetz aus.

$\oint_{C}\vec{H}\cdot \delta \vec{l}= I$ (1)

$\oint_{C}\vec{E}\cdot \delta \vec{l}=-\frac{\partial \Phi }{\partial t}$ (2)

Wo sind:

H [A/m] Intensität des magnetisches Feldes. (zeitlich veränderlicher Vektor, Richtung parallel zu Spulenachse).
I [A] freier Strom, der durch irgendeine Fläche an der Grenze C fließt. Gebundene Ströme werden vernachlässigt.
E [V/m] Intensität des elektrischen Feldes (zeitlich veränderlicher Vektor, Richtung tangential).
$\partial \Phi /\partial t$ [Wb/s] Veränderung der magnetischen Induktionsströmung in der Zeit

Fig. 1 – cylindrical workpiece undergoing induction heating, source [1]

Bild 1 – Zylindrisches Werkstück bei der Induktionserwärmung, Quelle [1]

Durch die Integration der Gleichung (1) um das auf dem Bild 1 markierte Rechteck herum und durch den Ersatz des Stroms durch das Produkt der Stromdichte J [A/m2] und der Fläche des markierten Rechteckes ergibt sich:

$(H_z+\delta H_z)l + H_z l=J_\phi \delta r l$ (3)

Die Stromdichte J in der tangentialen Richtung am Durchmesser r ist:

$J_\phi=-\frac{\partial H_z}{\partial r}$ (4)

Wir setzen das homogene Material mit der relativen Permeabilität $\mu_r$ und Resistivität $\rho [\Omega m]$ voraus. Für die Resistivität gilt das Ohmsche Gesetz in der Differenzialform:

$\vec{E}=\rho \vec{J}$ (5)

Für die magnetische Strömung durch die Fläche S, die senkrecht zu dieser Strömung ist, gilt:

$\Phi =\mu_0 \mu_r HS$ (6)

Durch die Integration der Gleichung (2) um den auf dem Bild 1 markierten Kreisring herum und durch das Einsetzungen der Proportionen (5) und (6) ergibt sich:

$2\pi (r+\delta r)\rho(J_\phi+\delta J_\phi)-2\pi r\rho J_\phi=-\mu_0 \mu_r 2 \pi r\frac{\partial H_z}{\partial t}\partial r$ (7)

Wir entfernen das unwesentlich kleine Glied $2\pi r\rho J_\phi$ und vereinfachen die Gleichung:

$\frac{\partial J_\phi}{\partial r}+\frac{J_\phi}{r}+\frac{\mu_0 \mu_r}{\rho} \frac{\partial H_z}{\partial t}=0$ (8)

Durch die Verbindung der Gleichungen (4) und (8) bekommen wir die 1D Gleichung des magnetischen Feldes der langen Spule:

$\frac{\partial^2 H_z}{\partial r^2}+\frac{1}{r} \frac{\partial H_z}{\partial r}-\frac{\mu_0 \mu_r}{\rho} \frac{\partial H_z}{\partial t}=0$ (9)

Die analytische Lösung dieser Gleichung steht in [1]. Ist die Verteilung der Intensivität des magnetischen Feldes im Werkstück gelöst, lässt sich mithilfe der Gleichung (4) die Stromdichte berechnen. Die Wärmeentwicklung im Material legen wir als die Joulesche Wärme fest. Die Hystereverluste lassen sich in bei meisten Anwendungsbereichen der Induktionserwärmung vernachlässigen. Die Verlustdichte p [W / m3] wird wie folgt berechnet:

$p=\rho J^2$ (10)

Folgende Grafen zeigen den Ablauf der gelösten Größen im erwärmenden Werkstück.


Bild 2 – Intensität des magnetischen Feldes	Bild 3 – Stromdichte	Bild 4 – Verlustdichte

Die Temperatur in dem Werkstück lösen wir mithilfe der 1D Gleichung für die nichtstationäre Wärmeleitung, die für die radiale Richtung diese Form hat:

$\frac{\partial T}{\partial t} = \frac{1}{r}\frac{\partial }{\partial r}\left ( \frac{\lambda r}{\rho c} \frac{\partial T}{\partial r} \right )+\frac{p}{\rho c}$ (11)

wo sind:
$T$ [°C] Temperatur
$t$ [s] Zeit
$\lambda$ [W/mK] Wärmeleitfähigkeit
$\rho$ [kg/m³] Dichte
$c$ [J/kgK] Spezifische Wärmekapazität
$p$ [W/m³] Verlustdichte

In die Berechnung werden die Verluste durch Ausstrahlung einbezogen als die Randbedingung für die Oberfläche. Die Verluste durch die Ausstrahlung entsprechen der vierten Potenz der thermodynamischen Temperatur der Oberfläche und bei der Erwärmung auf hohe Temperaturen können nicht vernachlässigt werden. Die Dichte der Wärmeströmung der Oberfläche wird mithilfe des Stefan-Boltmann-Gesetzes berechnet:

$j=\varepsilon \sigma \left ({T_{R}}^{4}-{T_{0}}^{4} \right )$ (12)

wo sind:
$j$ [W/m²] Dichte der Wärmeströmung durch die Oberfläche
$\varepsilon$ gegenseitiger Emissionsgrad des erwärmenden Körpers und der Umgebung (Ausmauerung)
$\sigma$ Stefan-Boltzmann- Konstante, ihr Wert ist 5,67×10^-8 Wm^-2K^-4 $T_{R}$ [K] thermodynamische Temperatur der Oberfläche des erwärmenden Körpers
$T_{0}$ [K] thermodynamische Temperatur der Oberfläche der Umgebung (Ausmauerung)

So kann die Verteilung der Temperatur im Werkstück am Ende der Erwärmung aussehen:

Fig. 7 – core-surface temperature profile

Bild 7 – Temperaturprofil Kern - Oberfläche

Referenzen:
[1] M. W. Kennedy, "Magnetic Fields and Induced Power in the Induction Heating of Aluminium Billets", 2013, ISBN 9789175018102
[2] E. Rapoport, Y. Pleshivtseva, "Optimal Control of Induction Heating Processes", 2006, ISBN 9780849337543