Preskočiť na obsah

Nestacionárna analýza vedenia tepla

Táto časť predstavuje časovú diskretizáciu a iteračnú metódu riešenia analýzy vedenia tepla v pevných telesách pomocou metódy konečných prvkov (FEM). Riadiace rovnice a okrajové podmienky kontinua nájdete v časti Rovnica vedenia tepla.

Diskretizovaná rovnica (východiskový bod)

Diskretizáciou rovnice vedenia tepla (rovnica vedenia tepla (gov_he_main)) Galerkinovou metódou dostaneme

\[\begin{equation} K T + M \frac{\partial T}{\partial t} = F \label{eq:2.4.8} \end{equation}\]

kde

\[\begin{equation} K = \int\left( k_x \frac{\partial N^T}{\partial x}\frac{\partial N}{\partial x} + k_y \frac{\partial N^T}{\partial y}\frac{\partial N}{\partial y} + k_z \frac{\partial N^T}{\partial z}\frac{\partial N}{\partial z} \right) dV + \int hc N^T N ds + \int hr N^T N ds \label{eq:2.4.9} \end{equation}\]
\[\begin{equation} M = \int \rho c N^T N dV \label{eq:2.4.10} \end{equation}\]
\[\begin{equation} F = \int Q N^T dV - \int q_s N^T dS + \int{hc} T c N^T dS + \int{hcTr} ({T+Tr}) ({T^2 + T r^2}) N^T dS \label{eq:2.4.11} \end{equation}\]
\[\begin{equation} N = (N^1, N^2, \ldots, Ni) \label{eq:2.4.12} \end{equation}\]

Tu \(K\), \(M\), \(F\) a \(N\) predstavujú maticu vedenia tepla (vrátane konvekčných a radiačných okrajových príspevkov), maticu hmotnosti, vektor tepelného zaťaženia a maticu tvarových funkcií. Definície symbolov materiálových vlastností (\(\rho\), \(c\), \(k_x, k_y, k_z\), \(Q\), \(hc\), \(hr\) atď.) sú uvedené v časti Rovnica vedenia tepla.

Časová diskretizácia a iteračné riešenie

Rovnica \(\eqref{eq:2.4.8}\) je nelineárna nestacionárna rovnica. Pri použití spätnej Eulerovej metódy na časovú diskretizáciu sa pri známej teplote v čase \(t=t_0\) vypočíta teplota v čase \(t=t_0+\Delta t\) pomocou nasledujúcej rovnice.

\[\begin{equation} K_{t=t_0+\Delta t} T_{t=t_0+\Delta t} + M_{t=t_0+\Delta t} \frac{T_{t=t_0+\Delta t} - T_{t=t_0}}{\Delta t} = F_{t=t_0+\Delta t} \label{eq:2.4.13} \end{equation}\]

Uvažujme zlepšenie vektora teploty \(T_{t=t_0+\Delta t}^{(i)}\), ktorý približne spĺňa rovnicu \(\eqref{eq:2.4.13}\), tak, aby sme získali presnejšie riešenie \(T_{t=t_0+\Delta t}^{(i)+1}\).

Na tento účel najprv vyjadríme vektor teploty nasledovne.

\[\begin{equation} T_{t=t_0+\Delta t}= T_{t=t_0+\Delta t}^{(i)} + \Delta T_{t=t_0+\Delta t}^{(i)} \label{eq:2.4.14} \end{equation}\]

Súčin matice vedenia tepla a vektora teploty, matica hmotnosti a súvisiace členy sa aproximujú nasledovne.

\[\begin{equation} K_{t=t_0+\Delta t} T_{t=t_0+\Delta t} = K_{t=t_0+\Delta t}^{(i)} T_{t=t_0+\Delta t}^{(i)} + \frac{\partial \big(K_{t=t_0+\Delta t}^{(i)} T_{t=t_0+\Delta t}^{(i)}\big) } {\partial T_{t=t_0+\Delta t}^{(i)} } \Delta T_{t=t_0+\Delta t}^{(i)} \label{eq:2.4.15} \end{equation}\]
\[\begin{equation} M_{t=t_0+\Delta t} = M_{t=t_0+\Delta t}^{(i)} + \frac{\partial M_{t=t_0+\Delta t}^{(i)}}{\partial T_{t=t_0+\Delta t}^{(i)}} \Delta T_{t=t_0+\Delta t}^{(i)} \label{eq:2.4.16} \end{equation}\]

Dosadením rovníc \(\eqref{eq:2.4.14}\), \(\eqref{eq:2.4.15}\) a \(\eqref{eq:2.4.16}\) do rovnice \(\eqref{eq:2.4.13}\) a zanedbaním členov druhého a vyšších rádov dostaneme nasledujúcu rovnicu.

\[\begin{equation} \bigg(\frac{M_{t=t_0+\Delta t}^{(i)}}{\Delta t} + \frac {\partial M_{t=t_0+\Delta t}^{(i)} } { \partial T_{t=t_0+\Delta t}^{(i)} } \frac{T_{t=t_0+\Delta t}^{(i)} - T_{t=t_0}}{\Delta t} + \frac{\partial \big(K_{t=t_0+\Delta t}^{(i)} T_{t=t_0+\Delta t}^{(i)}\big)} {\partial T_{t=t_0+\Delta t}^{(i)}} \bigg) \Delta T_{t=t_0+\Delta t}^{(i)} \\\ = F_{t=t_0+\Delta t} - M_{t=t_0+\Delta t}^{(i)} \frac{T_{t=t_0+\Delta t}^{(i)} - T_{t=t_0}}{\Delta t} - K_{t=t_0+\Delta t}^{(i)} T_{t=t_0+\Delta t}^{(i)} \label{eq:2.4.17} \end{equation}\]

Matica koeficientov na ľavej strane sa ďalej aproximuje pomocou nasledujúcej rovnice.

\[\begin{equation} K^{(i)} = \frac{M_{t=t_0+\Delta t}^{(i)}}{\Delta t} + \frac{\partial \big( K_{t=t_0+\Delta t}^{(i)} T_{t=t_0+\Delta t}^{(i)} \big)}{\partial T^{(i)}_{t=t_0+\Delta t}} = \frac{M_{t=t_0+\Delta t}^{(i)}}{\Delta t} + K_{T_{t=t_0+\Delta t}}^{(i)} \label{eq:2.4.18} \end{equation}\]

Tu \(K_{T_{t=t_0+\Delta t}}^{(i)}\) je dotyková matica tuhosti.

Napokon možno teplotu v čase \(t=t_0+\Delta t\) vypočítať iteráciami pomocou nasledujúcej rovnice.

\[\begin{equation} K^{(i)} \Delta T_{t=t_0+\Delta t}^{(i)} = F_{t=t_0+\Delta t} - M_{t=t_0+\Delta t}^{(i)} \frac{T_{t=t_0+\Delta t}^{(i)} - T_{t=t_0}}{\Delta t} - K^{(i)} T_{t=t_0+\Delta t}^{(i)} \label{eq:2.4.19} \end{equation}\]

Pri stacionárnej analýze sa konkrétne vykonávajú iterácie pomocou nasledujúcej rovnice.

\[ K_T^{(i)} \Delta T_{t=\infty}^{(i)} = F_{t=\infty} - K_T^{(i)} \Delta T_{t=\infty}^{(i)} \]
\[\begin{equation} T_{t=\infty}^{(i+1)} = T_{t=\infty}^{(i)} + \Delta{T}_{t=\infty}^{(i)} \label{eq:2.4.20} \end{equation}\]

Keďže sa pri nestacionárnej analýze na časovú diskretizáciu používa implicitná metóda, voľba časového prírastku \(\Delta t\) vo všeobecnosti nepodlieha obmedzeniu jeho veľkosti. Ak je však časový prírastok \(\Delta t\) príliš veľký, zvyšuje sa počet iterácií potrebných na konvergenciu. Vo všeobecnosti nadmerne veľký časový prírastok \(\Delta t\) zvyšuje počet iterácií. V implementácii automatické riadenie prírastku sleduje veľkosť vektora rezídua, znižuje \(\Delta t\) pri pomalej konvergencii a zvyšuje \(\Delta t\) pri malom počte iterácií (→ podrobnosti pozrite v časti Riadenie kroku).

Súvisiace témy

AI-assisted translation May contain errors Official docs Status