Page MenuHomec4science
Files F67269700

manual-heattransfermodel.tex
No OneTemporary
Actions

File Metadata

Created
Fri, Jun 21, 08:45

manual-heattransfermodel.tex
View Options

\chapter{Heat Transfer model}
The heat transfer model is a specific implementation of the \code{Model} interface
dedicated to handle the dynamic heat equation.
\section{Theory}
The strong form of the heat equation
can be expressed as
\begin{equation}
\rho c_v \dot{T} + \nabla \cdot \vec{\kappa} \nabla T = b
\end{equation}
with $T$ the temperature scalar field,
$c_v$ is the specific heat capacity,
$\rho$ is the mass density,
$\mat{\kappa}$ is the conductivity tensor and $b$ is the heat generation per unit of volume.
The weak form discretized with finite elements becomes
\begin{equation}
\forall i \quad
\sum_j \left( \int_\Omega \rho c_v N_j N_i d\Omega \right) \dot{T}_j
- \sum_j \left( \int_\Omega \vec{\kappa} \nabla N_j \nabla N_i d\Omega \right) T_j =
- \int_{\Gamma} N_i \vec{q} \cdot \vec{n} d\Gamma + \int_\Omega b N_i d\Omega
\end{equation}
with $i$ and $j$ being node indexes, $\vec{n}$ the field of normals of the surface
$\Gamma = \partial \Omega$.
To simplify, we can define the capacity and conductivity matrices
\begin{equation}
C_{ij} = \int_\Omega \rho c_v N_j N_i d\Omega \qquad \textrm{and} \qquad
K_{ij} = - \int_\Omega \vec{\kappa} \nabla N_j \nabla N_i d\Omega
\end{equation}
so that the discrete system to solve can be written
\begin{equation}
\mat{C} \cdot \vec{\dot{T}} = \vec{Q^{ext}} -\mat{K} \cdot \vec{T}~,
\end{equation}
with $Q^{ext}_i$ the consistent heat generated on node $i$.
\section{Using the heat transfer model}
At the present time an explicit dynamic analysis is available by
using the \code{HeatTransferModel} object.
It can simply be created like this:
\begin{cpp}
HeatTransferModel model(mesh, spatial_dimension);
\end{cpp}
while an existing mesh has been used (see \ref{sect:common:mesh}).
Then the model object can be initialized with:
\begin{cpp}
model.initFull("material.dat")
\end{cpp}
where a material file name has to be provided. This function will load the material
properties, and allocate/initialize the node and element vectors.
More precisely, the heat transfer model contains 4 \code{Vectors}:
\begin{description}
\item[temperature] contains the nodal temperatures $T$ (zero by default after the
initialization)
\item[temperature\_rate] contains the variation of temperature $\dot{T}$
(zero by default after the
initialization)
%\item[boundary] contains a \code{bool} value for each node
% specifying whether temperature is fixed or not. A Dirichlet boundary
% condition can be prescribed by setting the \textbf{boundary} value of a node
% to \code{true}. A Neumann boundary condition can only be applied
% if the \textbf{boundary} value of a node is \code{false}. If a
% node has a \textbf{boundary} value that is \code{false}, the
% \textbf{temperature}, \textbf{temperature\_rate} and
% \textbf{residual} are computed by the solve algorithm when relevant\todo{wait for till correction},
% otherwise these vectors contain the imposed values.
\item[boundary] contains a boolean value for each degree of freedom specifying
whether that degree is blocked or not. A Dirichlet boundary condition can be
prescribed by setting the \textbf{boundary} value of a degree of freedom to
\code{true}. A Neumann boundary condition can be applied by setting the
corresponding \textbf{boundary} value to \code{false}. The
\textbf{temperature} and the \textbf{temperature\_rate} are computed for all
degrees of freedom for which the \textbf{boundary} value is set to
\code{false}. For the remaining degrees of freedom, the imposed values (zero
by default after initialization) are kept.
\item[residual] contains the difference between external and internal heat generation. On
temperature imposed nodes, \textbf{residual} contains the support heat reactions. ($\vec{R} = \vec{Q^{ext}} -\mat{K} \cdot \vec{T}$)
\end{description}
Only a single material can be specified on the domain.
A material file allows to give the material properties in a text file (\eg material.dat) as follows:
\begin{cpp}
heat %\emph{name\_material}% [
capacity = %\emph{XXX}%
density = %\emph{XXX}%
conductivity = [%\emph{XXX}% ... %\emph{XXX}%]
]
\end{cpp}
where \code{capacity} and \code{density} are scalars, and the \code{conductivity} is specified as a
full tensor with $3\times 3$ components.
\subsection{Explicit dynamic}
The implemented explicit time integration scheme in \akantu uses a lumped capacity
matrix $\mat{C}$ (reducing the computational cost, see Chapter \ref{sect:smm}).
This matrix is assembled by
distributing the capacity of each element onto its nodes. The resulting $\mat{C}$ is
therefore a diagonal matrix stored in the \code{capacity} vector of the model.
\begin{cpp}
model.assembleCapacityLumped();
\end{cpp}
\index{HeatTransferModel!assembleCapacityLumped}
\note{At current time only explicit time integration with lumped capacity matrix
is implemented within \akantu.} \\
The explicit integration scheme is a \index{Forward Euler} {\it Forward Euler}
\cite{curnier92a}.
\begin{itemize}
\item Predictor: $\vec{T}_{n+1} = \vec{T}_{n} + \Delta t \dot{\vec{T}}_{n}$
\item Update residual: $\vec{R}_{n+1} = \left( \vec{Q^{ext}_{n+1}} - \vec{K}\vec{T}_{n+1} \right)$
\item Corrector : $\dot{\vec{T}}_{n+1} = \mat{C}^{-1} \vec{R}_{n+1}$
\end{itemize}
The explicit integration scheme is conditionally stable. The time step has to be
smaller than the stable time step which is obtained in \akantu as follows:
\begin{cpp}
time_step = model.getStableTimeStep();
\end{cpp}
\index{HeatTransferModel!StableTimeStep}
The stable time step is defined as:
\begin{equation}\label{eqn:smm:explicit:stabletime}
\Delta t_{\st{crit}} = 2 \Delta x^2 \frac{\rho c_v}{\mid\mid \mat{\kappa} \mid\mid^\infty}
\end{equation}
where $\Delta x$ is a characteristic length (\eg the inradius in the case of
linear triangle element), $\rho$ is the density, $\mat{\kappa}$ is the conductivity tensor
and $c_v$ is the specific heat capacity. It is
necessary to impose a time step that is smaller than the stable time step, for
instance, by multiplying the stable time step by a safety factor smaller than
one.
\begin{cpp}
const Real safety_time_factor = 0.1;
Real applied_time_step = time_step * safety_time_factor;
model.setTimeStep(applied_time_step);
\end{cpp}
%The initial temperature and temperature\_rate fields are, by default, equal to zero if
%not given specifically by the user (similar to the SolidMechanicsModel initial conditions
%\ref{sect:smm:initial_condition}).
The following loop allows, for each time step, to update the temperatures, residual and
temperature\_rate fields following the previously described integration scheme.
\begin{cpp}
for (UInt s = 1; (s-1)*applied_time_step < total_time; ++s) {
model.explicitPred();
model.updateResidual();
model.explicitCorr();
}
\end{cpp}
\index{HeatTransferModel!explicitPred}
\index{HeatTransferModel!explicitCorr}
\index{HeatTransferModel!updateResidual}
\index{HeatTransferModel!solveExplicitLumped}
An example of explicit dynamic heat propagation is presented in \\
\shellcode{\examplesdir/heat\_transfer/explicit\_heat\_transfer.cc}. \\
This example consists of a square 2D plate of \unit{1}{\squaremetre}
having an initial temperature of \unit{100}{\kelvin} everywhere but a none centered hot point
maintained at \unit{300}{\kelvin}. Figure \ref{fig:htm:explicit:dynamic} presents the
geometry of this case. The material used is a linear fictitious elastic material
with a density of \unit{8940}{\kilogrampercubicmetre}, a conductivity of
\unit{401}{\watt\per\metre\per\kelvin} and a specific heat capacity of \unit{385}{\joule\per\kelvin\per\kilogram}. The time step used is \unit{0.12}{\second}.
\begin{figure}[!htb]
\centering
\includegraphics[width=.4\textwidth]{figures/hot-point-1}
\hfill
\includegraphics[width=.4\textwidth]{figures/hot-point-2}
\caption{Initial temperature field (left) and after 15000 time steps = 30 minutes (right). The lines represent iso-surfaces.}
\label{fig:htm:explicit:dynamic}
\end{figure}

Event Timeline

c4science · Help