Thermische Simulation von Transport- und Lagerbehältern

Ein Transport- und Lagerbehälter für abgebrannte Brennelemente (T/L-Behälter) muss verschiedene Aufgaben und Anforderungen erfüllen. Die Hauptaufgaben eines solchen Behälters sind die sichere Abfuhr der Nachzerfallswärme, das Verhindern einer Kontamination der direkten Umgebung, das Sicherstellen einer konstanten Unterkritikalität sowie die Minimierung der ionisierenden Strahlung. Ein T/L-Behälter besteht im Wesentlichen aus einem Außenbehälter, einem Tragkorb zur Aufnahme der Brennelemente und einem Verschlusssystem, welches sich in der Regel aus zwei Deckeln zusammensetzt. In Abbildung 1 sind die Hauptkomponenten eines T/L-Behälters am Beispiel des CASTOR V/19 dargestellt.

Abbildung 1: Hauptkomponenten eines T/L-Behälters am Beispiel des CASTOR V/19 [1]

Als Resultat aus Fertigung und Nutzung der Behälter entsteht zwischen Tragkorb und Außenbehälter ein meist nur wenige Millimeter breiter, gasgefüllter Spalt, in welchem sowohl die Konvektion, die Wärmeleitung als auch die Wärmestrahlung am Wärmeübergang beteiligt sind. Dies führt dazu, dass die Modellierung und die thermische Simulation eines T/L-Behälters sehr zeit- und rechenressourcenintensiv sind. Weiterhin bedingt das geometrische Verhältnis von Spaltbreite zu Behältergröße eine sehr feine Vernetzung im Inneren des Spaltes, was den Aufwand zusätzlich erhöht.

Eine Möglichkeit, hier Abhilfe zu schaffen, ist das Umgehen der Vernetzung im Spalt durch den Einbau von Zwangsbedingungen. Diese Zwangsbedingungen müssen so formuliert sein, dass sie die Wärmeübertragung durch Strahlung, Konvektion und Leitung im Spalt hinreichend genau abbilden. Derartige Zwangsbedingungen wurden von Dinkel in [2] als „Thermische Spaltbedingung“ definiert. Diese integrieren analytische Wärmeübertragungsgleichungen in das Gesamtgleichungssystem der thermischen Finite-Elemente-Analyse (FEA). Die Gleichungen müssen dabei in einer Form vorliegen, die die Temperaturdifferenz zwischen den Knoten des Tragkorbs  T_i und den Knoten des Außenbehälters  T_a mit einer Konstante gleichsetzt [2]:

 T_i - T_a = K

Die Konstante 𝐾 wird dabei anhand der jeweiligen Wärmeübertragungsgleichungen bestimmt, sie darf lediglich keine Temperaturabhängigkeit aufweisen. Weiterhin muss für jede Form der Wärmeübertragung die aus der Behälterposition resultierende Spaltgeometrie berücksichtigt werden. Bei der Zwischenlagerung eines T/L-Behälters steht dieser aufrecht, was zu einer symmetrischen Spaltgeometrie führt. Für den Transport werden die Behälter jedoch auf die Seite gelegt, wodurch der Tragkorb nicht mehr konzentrisch im Außenbehälter liegt. Die verschiedenen Spaltgeometrien sind in Abbildung 2 dargestellt.

Abbildung 1: Zusammenhang der radialen Spaltgeometrie mit der Behälterposition

 

Wärmeleitung

Der durch Wärmeleitung übertragene Wärmestrom ergibt sich im Falle des stehenden Behälters zu [2, 3]:

 \dot{Q} = \frac{2*\pi*l*\lambda*(T_i -T_a)}{\ln{\frac{r_a}{r_i}}}

Dabei bezeichnet l die Behälterhöhe, \lambda die Wärmeleitfähigkeit des Fluids und r_a respektive r_i den Außen- und Innenradius des Spalts. Ein Umformen der Gleichung nach der Temperaturdifferenz ergibt die gesuchte Konstante [2]:

 K = \frac{\dot{Q}}{2*\pi*l*\lambda} * \ln{\frac{r_a}{r_i}}

Da diese Gleichung nur für einen Ringspalt gültig ist, muss der Wärmestrom bei der Transportkonfiguration jeweils zwischen zwei gegenüberliegenden Knoten separat berechnet werden. Hierfür müssen die Fläche  A der an der Wärmeübertragung beteiligten Oberfläche sowie die lokale Spaltweite  s_i berücksichtigt werden [3]. Mit Gleichung (1) folgt somit [2]:

 K_i = \frac{\dot{Q} s_i}{A*\lambda}

Konvektion

Im Falle der Konvektion, welche im gasgefüllten Spalt eines T/L-Behälters als freie Konvektion auftreten kann, spielt die Position des Behälters aus zwei Gründen eine wichtige Rolle: Zum einen ist, wie gehabt, die lokale Spaltweite zwischen zwei Knoten von Bedeutung, zum anderen muss hier ebenfalls die Richtung der wirkenden Schwerkraft berücksichtigt werden. Maßgebend für die thermische Spaltbedingung ist allerdings nicht die Kenntnis über das gesamte Strömungsfeld. Es reicht aus, zu wissen welchen Anteil die freie Konvektion an der Wärmeübertragung hat. Mit dieser Kenntnis kann die Wärmeleitfähigkeit des Fluids im Spalt derart modifiziert werden, sodass der Beitrag der Konvektion zur Wärmeübertragung in der Gesamtgleichung berücksichtigt wird. Dieses Vorgehen beruht auf der Annahme, dass die Wärmeströme von Wärmeleitung und Konvektion unabhängig voneinander betrachtet werden, was aufgrund des vergleichsweise geringen Anteils der Konvektion an der Gesamtwärmeübertragung eine zulässige Annahme ist. Die Nusselt-Zahl  N u als dimensionslose Kennzahl der freien Konvektion dient hier als Beurteilungskriterium [2]:

 \lambda_{LK} = Nu * \lambda_{fluid}

Die Wärmeleitfähigkeit unter Berücksichtigung von Wärmeleitung und Konvektion  \lambda_{LK} kann anschließend in Gleichung (3) bzw. (4) eingesetzt werden. Die Berechnung der Nusselt-Zahl erfolgt anhand verschiedener empirischer Korrelationen, welche für geometrisch einfache und technisch wichtige Fälle vorliegen. Im Fall des stehenden Behälters wird die Annahme getroffen, dass die Nusselt-Zahl über den Spaltumfang und die Spalthöhe hinweg konstant ist. Sie muss also nur einmal berechnet werden. Beim liegenden Behälter ist dies jedoch nicht mehr der Fall. Hier muss für jedes Knotenpaar entlang des Spaltumfangs, an welchem Konvektion zur Wärmeübertragung beiträgt, eine eigene Nusselt-Zahl berechnet werden. Für die genaue Herleitung der verwendeten Korrelationen wird an dieser Stelle auf die Arbeit von Dinkel [2] verwiesen. Im Fall des stehenden Behälters wird gemäß Literaturempfehlung der Maximalwert aus den folgenden Nusselt-Korrelationen verwendet [2]:

 Nu = \frac{C_1 * Ra * (\frac{l}{s})^2}{C_2 * (\frac{l}{r_a}^4*(\frac{r_i}{l}+[Ra * \frac{l}{s})^3]^{n_1}*(\frac{r_i}{l})^{n_2}}

Die Nusselt-Zahl ist also abhängig von der Rayleigh-Zahl  Ra . Die Bestimmung der Hilfsgrößen  C_1, C_2, n_1 und n_2 erfolgt tabellarisch in Abhängigkeit einer weiteren Hilfsgröße nach [2]. Um die asymmetrische Spaltgeometrie beim liegenden Behälter zu berücksichtigen, muss bei jedem Knotenpaar die lokale Spaltbreite  s_i beachtet werden. Die Nusselt-Zahl berechnet sich damit zu [2]:

 Nu =0,20*Ra^{0,25}*(\frac{r_1}{r_i})^{0,5} = 0,20*Ra^{0,25}*(\frac{r_i+s_i}{r_i})^{0,25}

 

Wärmestrahlung

Die Wärmeübertragung durch Strahlung zwischen zwei Oberflächen wird durch das Stefan-Boltzmann-Gesetz beschrieben [3, 4]:

 \dot{Q}_{12} = \sigma_{12} * A_1 * (T^4_1-T^4_2)

Die Indizes 1 und 2 repräsentieren hier die beiden betroffenen Oberflächen, wobei Fläche 1 die emittierende Fläche ist, im Falle des T/L-Behälters also die Mantelfläche des Tragkorbs. Die Strahlungskonstante  \sigma_{12} setzt sich den Emissivitäten  \varepsilon_1 und  \varepsilon_2 , den Oberflächen  A_1 und den  A_2 , dem Sichtfaktor  F_{12} sowie der Stefan-Boltzmann-Konstante \sigma = 5,67*10^{-8} W/m^2K zusammen [3,4].

 \sigma_{12} = \frac{1}{\frac{1-\varepsilon_1}{\varepsilon_1} + \frac{1}{F_{12}} + \frac{1-\varepsilon_2}{\varepsilon_2} * \frac{A_1}{A_2}}*\sigma

Der Sichtfaktor ist eine dimensionslose, rein geometrische Größe. Er bezeichnet den Anteil der Strahlung ausgehend von Fläche 1 der direkt auf Fläche 2 auftrifft [3, 4]:

F_{12} = \frac{1}{A_1} \iint\limits_{A_2} \left( \iint\limits_{A_1} \frac{\cos \psi_1 \cdot \cos \psi_2}{\pi \cdot r^2} \, dA_1 \right)dA_2

Hier sind \psi_1 und \psi_2 die Winkel zwischen den Normalenvektoren der jeweiligen Oberfläche zum Verbindungsvektor zwischen beiden Oberflächen. In Zwischenlagerkonfiguration muss der Sichtfaktor nur für einen einzigen emittierenden Knoten berechnet werden, da alle Knoten aufgrund der Symmetrie denselben Sichtfaktor aufweisen. Bei der Transportkonfiguration hingegen, muss der Sichtfaktor mehrfach berechnet werden, da die Symmetrieeigenschaft hier nur für die Hälfte des Spalts anwendbar ist.
In Gleichung (5) wird bereits die Wärmeleitfähigkeit definiert, welche sowohl die Wärmeleitung als auch die Konvektion für den Einbau als Spaltbedingung berücksichtigt. Damit der Einfluss der Wärmestrahlung zusätzlich betrachtet werden kann, wird ein Strahlungsleitkoeffizient  \lambda_S eingeführt. Dieser folgt durch Gleichsetzen der Wärmeströme aus (4) und (8), wobei die Wärmeleitfähigkeit als Strahlungsleitkoeffizient interpretiert wird [2]:

 \lambda_S = 4 * \sigma_{12} * s * (\frac{T_1+T_2}{2})^3

Die gesamte repräsentative Wärmeleitfähigkeit ergibt sich damit als Überlagerung der Wärmeleitfähigkeiten [2]:

 \lambda_{SLK} = \lambda_{LK} + \lambda_S

Eingesetzt in Gleichung (3) und (4) folgt damit die jeweilige Konstante für die thermische Spaltbedingung in Zwischenlager- bzw. Transportkonfiguration [2]:

 K= \frac{\dot{Q}}{2*\pi*l*\lambda_{SLK}} * \ln{\frac{r_a}{r_i}}

 K_i = \frac{\dot{Q} s_i}{A*\lambda_{SLK}}

 

Einbaumethode der thermischen Spaltbedingung

Für den Einbau von Zwangsbedingungen in der FEA existieren verschiedene Methoden, wobei die geläufigsten die Penalty-Methode, das Lagrange-Verfahren und das gestörte Lagrange-Verfahren sind. Beim Lagrange-Verfahren wird ein sogenannter Lagrange-Multiplikator 𝜆L verwendet, um Zwangsbedingungen in das Gesamtgleichungssystem der FEA einzubauen. Dadurch erhöht sich allerdings die Anzahl der Gleichungen, was sich wiederum negativ auf die Rechenzeit auswirkt. Vorteilhaft ist jedoch, dass die Wärmeleitfähigkeitsmatrix des Systems nicht manipuliert wird. Die Penalty-Methode basiert auf einer Multiplikation mit einem Strafparameter  \alpha_P [2,5]:

 (K + C^T * \alpha_P * C) * T = Q + C^T * \alpha_P * Q_P

Hierbei ist 𝑲 die Wärmeleitfähigkeitsmatrix des Systems, 𝑪 eine Koeffizientenmatrix mit aktiven Zwangsbedingungen, 𝑸 die Wärmemenge und 𝑻 die Temperatur.
Hierbei wird die Wärmeleitfähigkeitsmatrix manipuliert, die Anzahl an Gleichungen bleibt jedoch erhalten.
Das gestörte Lagrange-Verfahren implementiert beide Parameter, indem der Strafparameter  \alpha_{gL} verwendet wird um den Lagrange-Multiplikator  \lambda_{gL} zu manipulieren [2,5]:

 \left[  \begin{array}{cc}  K & C^T \\  C & -\frac{1}{\alpha_{gL}} \cdot I \\  \end{array}  \right] \cdot  \begin{bmatrix}  T \\  \lambda_{gL} \\  \end{bmatrix}  =  \begin{bmatrix}  Q \\  Q_\lambda \\  \end{bmatrix}.

Grundsätzlich sind alle genannten Möglichkeiten geeignet, um die thermische Spaltbedingung in eine Simulation einzubauen. Im Sinne der Zeitersparnis, was ein Grundgedanke der thermischen Spaltbedingung ist, ist die Verwendung der Penalty-Methode naheliegend.

Literatur

[1] https://www.gns.de/behaelter-equipment/brennelemente-haw/lwr/castor-v/castor-v19/
[2] C. Dinkel, “Integration analytischer Wärmeübertragungsberechnungsverfahren in das Finite-Elemente-System Z88 zur beschleunigten thermischen Bewertung von Transport- und Lagerbehältern für Brennelemente,” Dissertation, Shaker Verlag, 2019.
[3] Verein Deutscher Ingenieure, VDI-Gesellschaft Verfahrenstechnik und Chemieingenieurwesen (GVC), VDI-Wärmeatlas: Berechnungsblätter für den Wärmeübergang (VDI-Buch). Berlin, Heidelberg, s.l.: Springer Berlin Heidelberg, 2002.
[4] R. Marek and K. Nitsche, Praxis der Wärmeübertragung: Grundlagen – Anwendungen – Übungsaufgaben; mit 50 Tabellen, 50 vollständig durchgerechneten Beispielen sowie 126 Übungsaufgaben mit Lösungen; [auf CD: Lösungen der Übungsaufgaben und Programmbeispiele, 2nd ed. München: Fachbuchverl. Leipzig im Carl-Hanser-Verl., 2010.
[5] L. Nasdala, FEM-Formelsammlung Statik und Dynamik: Hintergrundinformationen, Tipps und Tricks, 3rd ed. Wiesbaden: Springer Vieweg, 2015.

 

Comments are closed.