These guys have never allreduced a timestep in a simulation where something can go aphysical.
I mean, sure, one can communicate aphysical conditions using out-of-band means. But it is so pretty when it just rides along for free as a NaN result.
Edit: I stand corrected having looked up authors. I still think having float effectively be Maybe[float] is a feature not a bug.
A closer analogue of NaN might be the x ("unknown") value in the three/four-valued logic of hardware simulation languages like Verilog. But the situation in Verilog, despite its challenges, is so much better since you don't have a hard discontinuity between "numbers" and "logical values". Numerical comparisons between unknowns can produce unknown logical values, and you have equations like true & unknown = unknown, true | unknown = true, ~unknown = unknown, etc, so you don't have the issue with NaNs where as soon as something crosses from the floating-point realm to the integer or logical realm, it has to commit to a definite value, with all the aforementioned issues.
All that said, NaNs were the right design choice given the constraints. But it's far from ideal.