After a little bit more thinking I changed the propagation function from
Heat(neighbor) += Heat(this) * decay / numberOfOpenNeighbors
to
Heat(neighbor) = Max(Heat(neighbor), Heat(this) - numOpenNeighbors * decay)
And now the decay factor is the valve for linear decreasing speed.
Intuitively you can see the ‘heat’ is decreasing linearly relative to ‘distance’ from heat source. But remember, we want a 1/distance^2 like curve for it, therefore we apply a post-processing step to the temperature we got for each grid. That is:
(1) Distance(grid) = Heat(source) - Heat(grid) (2) Factor(grid) = 1 / (Distance(grid) * Distance(grid) + 1) (3) Factor(grid) = Factor(grid) - 1 / (Heat(source) * Heat(source) + 1) (4) Heat(grid) = Heat(source) * Factor(grid)
By applying step (2), we generated a 1/distance^2 like function. (The +1 term on the denominator is to avoid divide by zero at source point.) But you can see that even a grid originally with Heat = 0 will be assigned with a positive factor. Step (3) is responsible for eliminating those residues.
Here I present the test results:


Compared with the result of last post, the result of this new mathematical model yields more reasonable area of influence, as well as easier-to-tune heatmap parameters.
Hope it helps you somehow!
PS: this algorithm takes 1~2 ms to run on my Unity3D scene, with grid map size ~60*40, which is acceptable for my current need, but not as good as expected. Later post I might talk about my experiments on optimizing the algorithm.