Lazy evaluation of a continuous field of signals #600
Conversation
I'm not really pleased with the fact I had to disable the `display_tile_overlay` system when using the signal overlay but this can be cleanup later. Color scheme has been updated to look better with the new terrain and shadows, and the overlay material has been made transparent a bit.
This only takes into account sources for now. This means that it only reaches part of the goal set in Leafwing-Studios#511: 1. It only tracks signal emitters. 2. It does determine an equation for each emitter. The equation simply is the solution of the diffusion equation for that emitter: `1/(1+D*t) * 1/sqrt(4*pi*k*t) * exp(d^2 / (4*k*t))` where: - `D` is a decay rate constant, and `1/(1+D*t)` is a decay term. - `k` is a diffusivity constant. - `d` is the distance between the emitter and a tile. - `t` is the time elapsed since the signal was emitted. 3. Contributions of each emitter are then integrated to compute the actual field value. This works because of the linearity of this subset of the problem, but may not hold when tackling e.g. barriers or signal modifiers (I haven't yet investigated it). 4. It lazyly looks up the field value when needed, but I suspect this is pretty costly when displaying the signal overlay. I don't know if and what we'd want to do about it yet.
Our application is in 2D. The 1D fundamental solution was used previously.
|
I did not have the time to push this draft as soon as I wanted, but here it is. There's 2 issues right now that are on top of my priority list:
Then next up probably is a bit of clean up and design talk, before I attempt tackling obstructions. |
alice-i-cecile
added
performance
|
Very well done! ⭐ I love that you implemented the debug overlay aswell! ✨ Would you be interested in a review? |
Yes please ! There's a lot to improve, so I'd gladly get feedback early to go in the right direction 😄 |
|
Super cool! I'm going to start by splitting out the signal-visualization work: this is valuable regardless of the approach used and I think we can handle it more robustly in its own PR :) |
Co-authored-by: Austreelis <[email protected]>
Co-authored-by: Austreelis <[email protected]>
* Extract initial signal viz code from #600 Co-authored-by: Austreelis <[email protected]> * Move tile selection code into infoviz module * Store an Option for visualized signal type * Use a single resource to store overlay info * Refactor material fetching into a method * Log error, don't panic * Reparameterize color ramp scaling for clarity * Fix math * Tweak max signal strength for better visualization * Use N_COLORS consistently * Display color ramp in legend * Organize color palette constants * Don't screw up the base color * Get tile interactions and signal viz to play nice together * Create color ramp via HSLA color interpolation * Store the legend in the right orientation * Tweak overlay colors for now * Fix math for normalizing signal strength * Attempt to do UI layout * Display currently visualized signal * Cycle signal types randomly * Clippy --------- Co-authored-by: Austreelis <[email protected]>
| use super::{SignalStrength, DECAY_RATE, DECAY_THRESHOLD, DIFFUSIVITY}; | ||
|
|
||
| #[derive(Clone, Copy, Debug)] | ||
| pub struct SignalEmission { |
There was a problem hiding this comment.
Doc strings for the type as a whole please :)
| /// | ||
| /// where: | ||
| /// - ϕ is the density of the diffusing material. In emergence's case, ϕ is the signal strength. | ||
| /// - D is the diffusivity of the environment. For emergence, this is a constant for now, making |
There was a problem hiding this comment.
| /// - D is the diffusivity of the environment. For emergence, this is a constant for now, making | |
| /// - D is the diffusivity of the environment. For Emergence, this is a constant for now, making |
| /// | ||
| /// The fundamental solution solves a problem with no boundary conditions and | ||
| /// an initial density given by the [Dirac delta function](https://en.wikipedia.org/wiki/Dirac_delta_function). | ||
| /// The general solution, for any initial condition is a convolution with the fundamental solution: |
There was a problem hiding this comment.
| /// The general solution, for any initial condition is a convolution with the fundamental solution: | |
| /// The general solution for any initial condition is a convolution with the fundamental solution: |
| /// - y is a surface element in the same space as x. | ||
| /// | ||
| /// Practically, we don't want to deal with costly convolutions, so we'll try to avoid them. | ||
| /// Because signals are emitted in pulses, and convolution is distributive, we can separate each |
There was a problem hiding this comment.
Do we want signals to be implemented in pulses? I think ideally these would be better modeled as plateaus, to avoid situations where wave fronts propagate away from the actual source.
| /// s(x, t) = ∫ Φ(x−y, t-T0) e0(y) dy + ∫ Φ(x−y, t-T1) e1(y) dy + ... + ∫ Φ(x−y, t-Tn) en(y) dy | ||
| /// | ||
| /// where: | ||
| /// - en is the initial condition of the nth emission. |
There was a problem hiding this comment.
I think that in comments, this is clearer as e_n
| /// s(x, t) = Φ(x-y0, t-T0) + Φ(x-y1, t-T1) + ... + Φ(x-yn, t-Tn) | ||
| /// | ||
| /// where yn is the location of the nth signal emission. | ||
| // NOTE: One equation per signal type implies they don't interact together. |
There was a problem hiding this comment.
This should be a doc comment on emissions
| }); | ||
| } | ||
|
|
||
| /// Sets the time at which the equation will be evaluated to the given `current_time`. This |
There was a problem hiding this comment.
Nit: the explanatory text belongs in a new paragraph; this plays nicer with IDE tooltip hovers.
| // negligible being more for fast decaying signals (we'll trim it earlier), and | ||
| // less for slowly decaying ones (we'll trim it later) | ||
| Some(SignalEmission { time, strength, .. }) | ||
| if DECAY_THRESHOLD * (1.0 + DECAY_RATE * (self.time - time)) > strength => {} |
There was a problem hiding this comment.
I'm having trouble following the math here when skimming; I think intermediate variables will make this clearer.
|
Very well-made, and fascinating math. Main conclusions:
|
|
About 4: I'm guessing this is done by ants to follow gradients of signal, right ? If so caching results for the time of one signal pass can probably be an easy win, given how ants tends to be close together. |
|
About 4: I'm guessing this is done by ants to follow gradients of signal, right ? If so caching results for the time of one signal pass can probably be an easy win, given how ants tends to be close together. Correct :) |
alice-i-cecile
added
the
controversial
|
I'm unlikely to have time to work on this PR before after the bevy jam, but I'm still looking at it. If anyone has suggestions/changes to add, I can probably take a quick look. And I'll continue after the Jam :) |
Fixes #511
These patches implement an algorithm solving the diffusion equation "on-demand". It is still only a draft, and needs optimization to perform better than the fixed-time step diffusion passes.
Closed form of the solution to diffusion
As pointed out by @RobWalt here, solving the simplified diffusion equation in our case boils down to computing a convolution of a so-called fundamental solution and initial conditions. The fundamental solution is the solution with no border condition and a Dirac delta as initial conditions.
The approach I took may be different than what was outlined in the comment linked above, so I'll develop a bit further. I'll ignore border conditions (obstruction by structures) for now.
It would make sense to model emission of signals as Dirac delta functions (at the instant they are emitted), which then get diffused. This means that solving for 1 emission is trivial, this is just the fundamental solution:
Which we can adjust by keeping track of when and where signals are emitted:
u(x, t) = ϕ(x-x0, t-T0)For several emissions (potentially emitted at different times), we could solve for one emission while the second hasn't been emitted, then do the complete convolution, and so forth for the following ones. However, convolutions are distributive, so the solution is the sum of the contributions of each emissions, no matter when they were emitted.
So generally for
nemissions:u(x, t) = ϕ(x-x0, t-T0) + ϕ(x-x1, t-T1) + ... + ϕ(x-xn, t-Tn)This allows us to model the field exactly with 2 parameters per signal emissions.
Expressing decay
The fundamental solution of the simplified diffusion equation we use represents a "density" of the diffused material. We could see it as pheromone concentration. The field being a density means we can integrate over its surface to get a total "amount" of pheromones. This "volume" is constant for the fundamental solution previously described, which means no signal strength is lost, it only gets displaced.
We could interpret decay as a gradual loss of signal strength, which would need us to use a fundamental solution that keeps the same shape as without decay, but with a decreasing volume. This is actually quite easy, we just need to add a dividing factor to the fundamental solution:
The algorithm
We use 3 types:
SignalEmission, representing a single emission pulse: its location, the time when it was emitted, and its strength.SignalEmissions, and the current time.Signals, mapping signal types to their equation.Every time the
update_signalssystem runs, we increment the current time of eachDiffusionEquation. We then pruneSignalEmissions which are considered negligible: any emission whose strength is smaller thanTh (1 + D t)(whereDis the decay rate,tthe time since emission, andThsome threshold) is removed from the equation. This is equivalent to pruning emissions whose volume under the curve is less thanTh.To emit a signal, we simply add a
SignalEmissionto theDiffusionEquation, with emission time set to the current time.When a signal is accessed, we get the associated
DiffusionEquationand iterate over its emissions. We compute the fundamental solution for each emission and accumulate them before returning the result. There is one subtlety: The distance is computed in tiles (using hexx'Hex::distance_tofunction), meaning that contour lines of equal signal strength are hexagons instead of circles.Possible improvements
Some of the features the previous signals implementation had, and some that may be needed/nice later, are not handled yet by this algorithm:
The amount of signal emissions stored seems to be the bottleneck, making this slower than the previous signal diffusion. But I'm pretty sure we could prune way more aggressively emissions, and that past some time their shape is so similar we can approximate a bunch with just one.
There's also some other things to investigate:
Todos
In no particular order:
pubvspub(crate)vs private)