Thermal center of a set The current Wikipedia article on the heat equation includes a gif which shows how the solution to the heat equation evolves over time if the initial data is the indicator function of a domain.

Suppose $w(t,x)$ is the solution to the heat equation $w_t=\Delta w$ satisfying $w(0,x)=\mathbb{1}_A(x)$ for a nice set $A$. Suppose further that the function $w(t,x)$ (as a function of $x$) attains a unique maximum at $x=a(t)$ for each $t$ in some interval $(0,\varepsilon)$, and we can define the thermal center of $A$ to be the limit $a=\lim\limits_{t\to0^+} a(t)$. (I made up this concept.)
Can this be given an integral formula, or integral characterization? Or some other nice and interesting characterization? For example, maybe it happens to be the center of mass of $A$, or something. Or, what is it for, say, a triangle?
 A: TL;DR In any dimension the position(s) of the thermal center(s) belong to the set
$$\mathbf{x}_{th}=\text{argmax}_x \min_{r\in \partial A} d(x,r)$$
Here is a sketch of a proof in $d=2$. (To extend the idea to arbitrary $d$ one just needs to replace $C_r^x\to (S^{d-1})_r^x$ and make appropriate adjustments to the measures and estimates, but I don't think this precludes the proof from going through).  Assume that the set $A$ is compact. Define the density of neighbors of a point $x$ at distance $r$ from that point to be
$$\rho(r;x)=\frac{dN}{dr}(r;x)=r\mu(C^x_r\cap A)$$
where $C_r^x $ is a circle of radius $r$ centered at $x$ and $\mu(C_r^x)=2\pi$. This function in a certain sense measures how many neighbors the point $x$ has a certain distance away. Note that by definition
$$\int_0^{\max _{r\in A}d(x,r)}\rho(r;x)dr=\mu(A)$$
EXAMPLE: For a circle of radius $R$ centered at the origin, the density of neighbors can be computed to be
$$\rho(r;x)=\begin{Bmatrix}2\pi r&r\in (0,R-|x|)\\ 2r\arccos\left(\frac{|x|^2+r^2-R^2}{2|x|r}\right)&~~~~~~~~~r\in (R-|x|,R+|x|)\end{Bmatrix}$$
Critically, notice here that the neighbor density for any two points in the interior of $A$ is the same up until some radius $L=\min\{\min_{r\in A}d(x_1,A), \min_{r\in A}d(x_2,A)\}$. The reason why this is important is because the idea is that in the limit $t\to 0^+$ the difference in value between two points is determined by the difference in the local environment of the two points. In layman words, the furthest away you are from the boundary, the more points you have around you and therefore you get a bigger contribution in the limit.
With this notation established, notice now that the solution to the initial value problem posed in the question is given by
$$w(x,t)=\frac{1}{4\pi t}\int_{0}^{\max _{r\in A}d(x,r)}e^{-\delta^2/4t}\rho(\delta;x)d\delta$$
Finally, define $r_0=\max_x\min_{r\in \partial A}d(x,r)$. For an arbitrary point $x$, define $r_x=\min_{r\in \partial A}d(x,r)$ and $R_x=\max\{\max_{r\in \partial A}d(x,r), \max_{r\in \partial A}d(\mathbf{x}_{th},r)\}$. It should be clear that $r_x\leq r_0$. Consider the difference
$$4\pi t(w(\mathbf{x}_{th},t)-w(x,t))=\left(\int_{r_x}^{r_0}+ \int_{r_0}^{R_x}\right)e^{-\delta^2/4t}(\rho(\delta;\mathbf{x}_{th})-\rho(\delta;x))d\delta:= I_1+I_2$$
One can show easily that $I_1\geq 0$ and hence establish the bounds
$$C(x) e^{-r_0^2/4t}\leq I_1 \leq 4\pi t(e^{-r_0^2/4t}- e^{-r_x^2/4t})\\|I_2|\leq 4\pi t (e^{-r_0^2/4t}-e^{-R_x^2/4t})$$
for $C(x):=\int_{r_x}^{r_0}(2\pi \delta- \rho(\delta;x))d\delta$ independent of $t$. Thus, in the limit $t\to 0^+$ it can be shown that $|I_2|/I_1\to 0$, and it can also be shown that for a small enough time
$$ t < A\epsilon C(x)~,~ \epsilon<\min\left(1, \frac{W^2-r_0^2}{4}\right)$$
where $A$ a positive constant, and $W$ denotes the radius of the smallest circle containing $A$.  Hence, it must be true that $I_1+I_2>(1-\epsilon)I_1>0$. It follows immediately that
$$\forall t<\epsilon(x):= AC(x)\epsilon~~,~~ w(x,t)\leq w(\mathbf{x}_{th},t)$$
By definition, it is true that the function $C(x)$ is always bounded and non-negative on compact sets and zero only when $x\in \mathbf{x}_{th}$. If we restrict $x$ to only consider $x\notin \mathbf{x}_{th}$, by the estimate above it is seen that for any of those $x$, there exists a finite time interval for which the last statement holds, which allows us to conclude that the thermal set defined above contains the maximum in the limit $t\to 0$.
NOTE: Obviously when $\mathbf{x}_{th}$ contains only one element, then that is the unique thermal center. When it doesn't, the thermal center is determined by subleading contributions to the integral above. (My intuition is that thermal centers can only be non-unique when the set possesses a symmetry group that allows for the existence of multiple members of the thermal set with the same environment. In any other case, there must be a unique thermal center.)
COMMENT 1: You can think about the last step in the proof as follows. For any small enough finite non-zero time $t=T$ it guarantees that there is a set of points in the set $x\notin \mathbf{x}_{th}$ for which $w(x,T)<w(\mathbf{x}_{th},T)$ (the points for which the condition $C(x)>T/A\epsilon$ is fulfilled). As $T\to 0$, the set shrinks until it becomes the empty set.
COMMENT 2: There exist sets $A$ where $\mathbf{x}_{th}$ can be relatively large. An example is given by the set $$A=C_R^{(-L,0)}\cup ([-L,L]\times[-R,R])\cup C_R^{(L,0)}$$
By construction it admits the thermal set
$$\mathbf{x}_{th}=[-L,L]\times\{0\}$$
However, it appears that the thermal center is the origin, which can be justified by inspecting $\rho(r;x\in \mathbf{x}_{th})$.
