Add smooth_clamp function

[SouthEndMusic]

Fixes #1918.

The function I implemented is somewhat different than the one suggested in the issue. In the clamp function I replaced the sections on [lo - d, lo + d] and [hi - d, hi + d] C1 smooth with a quadratic polynomial section, where d = a (hi - lo) / 2, so 0 <= a <= 1 is the relative smoothing parameter. The only edge case this cannot handle smoothly is when isinf(lo) xor isinf(hi).
I haven't noticed any significant performance improvements with this but it's quite little regret.

[Huite]

This looks like a nice implementation! I'd replace α by δ, and omit the default value.
Maybe switch the conditionals around so it matches the behavior of the regular clamp function when lo and hi are provided "wrongly".

[Huite]

If either hi or lo are inf, we just get the hard clamp. Ideally though, I think you'd want the clamping one on side not to be influenced on what happens on the other side, i.e.:

smooth_clamp(x_near_lo, lo, inf) == smooth_clamp(x_near_lo, lo, hi)

Thinking about this, why α instead of δ? δ works regardless of finite lo or hi. Arguably, δ is more intuitive as well, since your docstring explicitly starts off with δ too:

Function that goes smoothly from lo to hi in the interval [lo - δ, hi + δ]

And then defines it in terms of α. So in use, I have to do the α arithmetic first.

I think in terms of general use, we'd mostly have an absolute range in mind? If a relative zone is useful, it can be either provided in the call; alternatively you could support both arguments in an abstol, reltol kind of way.

Having α makes it easier to choose a default value, but... we probably don't want a default value here! The width of the clamping can make a big difference. I.e. in this line:

q = smooth_clamp(q_unlimited, -max_flow_rate[node_id.idx], max_flow_rate[node_id.idx])

The smoothing width should probably be explicit.

[Huite]

Do we want to check for lo > hi, or 0 <= α <= 1?

...I just tested, the clamp function doesn't care about lo and hi, so it's probably okay to leave unchecked, see also its source:

function clamp(x::X, lo::L, hi::H) where {X,L,H}
    T = promote_type(X, L, H)
    return (x > hi) ? convert(T, hi) : (x < lo) ? convert(T, lo) : convert(T, x)
end

The behavior isn't the same though if you fill in values "wrongly":

julia> smooth_clamp(0.01, 1.0, 0.0,α=0.01)
1.0
julia> clamp(0.01, 1.0, 0.0)
0.0

Might be best in some sense to have it mimick the lo vs hi behavior.

[visr]

Putting this back in draft for now. @SouthEndMusic indicated he didn't measure a performance gain from this. Though I can also imagine this only hits particular models, or that there is more work than just this discontinuity needed to make the whole thing smooth.