Heat equation (steady state) with boundary conditions at domain edges and inside the domain

Question

I'm trying to solve the steady state of a heat equation problem in 2D $\Delta u = 0$ (3D also), with the method of solving the huge system of equations that arises from the discretization of the domain. My boundary conditions are different from any example that I've been able to find.

If my domain is a 2D square, I want to put a fixed temperature (say $T_e$) all around the square edges. However, I want to force a point in the center of the square to be at a different temp from the edges, $T_c, T_c > T_e$. Think of a square metal plate being heated at the center by a soldering iron. I want to calculate the steady state heat in such a scenario.

To the best of my knowledge it seems that this problem is not really a heat equation problem: having a point at the center of the domain with a local maxima temperature goes against the properties of the heat equation...

I'm obtaining results, but it seems that heat does not propagate linearly from the center to the edges, but rapidly decaying. The problem is the same if moved into the 3D domain (a sphere/cube that is heated from it's center to its exterior edges, and my results are the same.

Maybe I'm using a wrong model? or a wrong tool to solve this steady state problem? Any hints or directions will be highly appreciated.

Best regards

Answer 1

You need to add the condition that $u(c) = T_c$ as an explicit constraint via a Lagrange multiplier. For simplicity, assume that $T_e=0$ (otherwise construct a $u_e$ such that $u_e|_{\partial\Omega} = T_e$ and add it to the solution with homogeneous boundary conditions, shifting the constraint by $T_c-u_e(c)$).

You know that the solution of $-\Delta u = 0$ solves the minimization problem $$ \min_{u\in H^1_0(\Omega)} \int_\Omega |\nabla u|^2 \,dx. $$ To enforce the constraint $u(c)=T_c$, you consider the Lagrange functional $$ \max_{\lambda\in\mathbb{R}} \min_{u\in H^1_0(\Omega)} \int_\Omega |\nabla u|^2 \,dx + \lambda(u(c)-T_c). $$ The necessary optimality conditions for this problem are obtained (formally) by differentiating with respect to $u$ and $\lambda$: \begin{align} \int_\Omega \nabla u\cdot \nabla v + \lambda v(c)&= 0 \quad\text{for all } v\in H^1_0(\Omega)\\ (u(c) - T_c)\mu &=0 \quad\text{for all }\mu\in\mathbb{R}. \end{align} (Note here that except in 1D, functions in $H^1(\Omega)$ are not continuous and thus do not have point values. You probably should replace the point value with an average value over a small subdomain. For a finite element discretization, you move to a subspace of piecewise polynomials, where point values make sense.)

For details on saddle point problems for partial differential equations, see Dietrich Braess, Finite Elements, Cambridge University Press, 2007, Chapter III.4.

Answer 2

From a mathematical viewpoint you can't do what you are trying to do: the problem $$ \Delta u = 0 \qquad \text{in } \Omega, \\ u = T_e \qquad \text{on } \partial\Omega, \\ u(c) = T_c \qquad \qquad\qquad\text{} $$ does not, in general, have a solution. This is because, mathematically speaking, the space of solutions of the Laplace equation, $H^1$, consists of functions that may not be continuous and so it doesn't make sense to take point values such as $u(c)$.

An alternative explanation of this problem would be to realize that if $T_c>T_e$ then there will be heat flowing from $c$ to the boundary, and consequently someone will have to put that heat into the body at $c$. But putting a finite amount of heat into a single point (of area zero) means to have a device that can produce an infinite heat flux -- hardly realistic.

The solution to what you are trying to do is to assume that the area at which you are able to impose a temperature $T_c$ is finite, i.e., to consider a domain that has a hole of finite size at the center.

Answer 3

Based on some of the previous answers, I think you are faced with a few choices in your modeling. If you are forced to use a cartesian unskewed grid (not clear from your question what discretization you are using), it's going to be a little tougher because as said previously, you need to cut off a finite area of a few cells in the center and you might get corner issues once you try to solve this.

If you can do arbitrary grids, may I suggest a circular one (not a square) which means you can use the symmetry properties to reduce the dimensionality of your problem. You then can compare your results with an analytical solution of the form:

$T(r) = (T_1 - T_0)\frac{\ln{\frac{r}{r_0}}}{\ln{\frac{r_1}{r_0}}}+ T_0$

And you can see that as $r_0 \rightarrow 0$, the temperature at the origin becomes undefined because your heat flux (positive or negative) would become infinite.

Answer 4

I am not sure about interpreting $\nabla u=0$ as heat spreading, given that it corresponds to the steady state of the time dependent (parabolic) equation $u_t=\nabla u$. For me the meaning is rather closer to an interpolation of the values specified at the boundary of the domain of definition $\partial \Omega$. These values definitely determine the function $u$, so if they are set to zero on a ball of radius 1, the function is identically zero. This is consistent with the formula of the solution to $\nabla u=0$ with $u=f$ at the boundary, given by:

$u(x)=\frac{-1}{\omega_n} \int_{|s|=1} \frac{-1+|x|^2}{|s-x|^n} f(s) ds, x \in \Omega$

for $| \cdot |$ the norm of the Euclidean space.

Moreover, having an interior point with maximum value is against the maximum principle holding for all armonic functions.

Thus, I am not sure if a boundary problem on a unit ball with constant temperature at its boundary and another temperature at an interior point it is well formulated for $\nabla u=0$. At least the interior point should be removed, so the domain becomes a sort of spherical band.

In other words, given that the solution to $\nabla u=0$ depends on the function values set on the boundary of the domain and on the boundary shape itself, what are these in the case of the fundamental solution (given by a radial solution)?. Further, is there any categorization of generic (fundamental) solutions depending on boundary shapes and values?