I am trying to model the heat equation with heat source and Robin boundary conditions, i.e. the system:
\begin{align}
T_t\;=&\;\alpha\Delta T+\frac{1}{c_p\rho \text{Vol}(\Gamma)}1_{\Gamma}(\boldsymbol{x})w(t). && \boldsymbol{x}\in \Omega, t\in I\nonumber \\[1mm]
\nabla T \cdot \overset{\to}{n}_{\boldsymbol{x}}=&-K(\boldsymbol{x})(T(\boldsymbol{x},t)-T_{\text{outside}}(t)) &&\boldsymbol{x} \in \partial \Omega, t\in I \label{eq:boundary_condition}\\[1mm]
T(\boldsymbol{x},0)\;=&\;T_0(\boldsymbol{x}) && \boldsymbol{x} \in \Omega
\end{align}
in $3$ dimensions using the finite difference method. My question is the following:
The ovens (heat source) becomes unreasonably hot. Any ideas why?
Here are some more details. Our domain $\Omega$ is a rectangular prism representing a room, and the heat source $\Gamma$ is two electrical ovens modeled as two (disjoint) rectangular prisms in the interior of $\Omega$. We have discretized the room into nodes $p_{i,j,k}$, and we imagine that each node sits in the centre of a cube of side lenght $h$. In our equation, $w(t)$ is the oven wattage- a given function, and $1_\Gamma$ is the indicator function of the ovens. Our explicit FTCS-scheme is given from
\begin{align*}
\frac{T_{i,j,k}^{m+1}-T_{i,j,k}^{m}}{\Delta t}\doteq&\; \alpha \frac{T_{i+1,j,k}^m-2T_{i,j,k}^m+T_{i-1,j,k}^m}{h^2}+\alpha \frac{T_{i,j+1,k}^m-2T_{i,j,k}^m+T_{i,j-1,k}^m}{h^2}\\[1mm]
+\;& \alpha \frac{T_{i,j,k+1}^m-2T_{i,j,k}^m+T_{i,j,k-1}^m}{h^2}+\frac{1}{c_p\rho \text{Vol}(\Gamma)}1_{\Gamma}(p_{i,j,k})w(t_m).
\end{align*}
Example:
If we run only a few $(1-20)$ time-steps, the temperature in the oven nodes remain reasonable, but not if we run over a longer period. For example, if the parameters are
\begin{align}
\text{room dimensions}\;=&\;4m\times 6m\times 2m && \text{time length}=120\text{sec}\\
h=\text{spatial step size}\;=&\;1/15 m && \Delta t=1/100\text{sec}\\
w(t)\;\equiv&\; 750\text{watts} && \alpha=1.9\cdot 10^{-5}m^2/s\\
T_{\text{outside}}(t)\;\equiv&\;253.15\text{Kelvin}=-20\text{Celsius}
\end{align}
Then the temperature in the oven nodes end up at $1408$Kelvin=$1135$celsius after these two minutes.
What I have tried
-
Making the mesh finer, so that the ovens are ensured to have at least $3-5$ nodes of thickness.
-
Giving the ovens larger volume.
-
Changing the ovens, so that they have larger surface area
-
Moving the ovens further away from the boundaries of the room, to avoid the problem being caused by the boundary conditions on the walls.
-
Turning down the oven wattage.
-
Using the corresponding equation one gets from temperature-dependendt density, heat capacity, diffusivity and conductivity of air.
Unfortunately, I still get the same problem in all cases (also if i change the step sizes, time step size and the other parameters.)
Some thoughts
Some tests have shown that: The ovens do indeed output the correct amount of energy into the room. The boundary conditions on the room also work as expected. If the initial temperature equals the outside temperature (both constant), and the ovens are turned off, then the temperature in the room remains constant as expected.
Our Courant/ Mesh Fourier number $\frac{\alpha \Delta t}{h^2}$ is quite small, so I don't believe that the step-sizes is the problem. (Although I suspected this number was in fact too small. Making it larger did not improve the situation).
Two other possible causes may be that the heat is not able to diffuse away from the oven fast enough, or that the problem is somehow caused by the heat-term. As mentioned, I tried making the diffusivity temperature-dependent, but it did not work. Also, it seems unreasonable that this adjustment should be necessary, since the heat is evenly distributed over the oven, since many panel ovens are mostly air anyway, and since some papers argue that realistic modelling can be obtained even with diffusivity assumed constant (e.g. this paper.)
It should be mentioned though that the result may be close to what one would expect if no heat leaves the ovens at all.
Best Answer
The results were indeed reasonable. The problem is that diffusivity of heat in air is a slow process, and that we should have included advection in our model. This can be done e.g. by modelling the advection-diffusion equation. Due to time constraints, we chose to go with the simpler approach of distributing the heat produced by the ovens uniformly over the room, i.e. we put $$\Gamma:=\Omega.$$ This gave much more reasonable results. We also lowered the wattage, as it was very high compared to the size of the room.
One could question how realistic it is to make this simplification. We think this assumption is realistic enough when the room in question has a reasonable size. This is because the heat from the ovens will quickly flow upwards around the ovens and create so-called "Bénard cells" in the room. It is also possible that the electrical panel ovens have fans built into them, which also makes the heat distribute faster.
Those who are doing similar projects, with similar struggles, can see how we implemented our method in Python, as well as our results, here.