Best way to forcibly passify a result using ForwardDiff?

Specific Domains
, ,
1

I’m working on the implementation of a time-stepping method and would like to force that part of the computation always be a passive value (in the AD sense).

The function in question is a simple broadcasted computation, followed by a weighted average:

const _WENO_OPTIMAL_STENCIL_WEIGHTS = @SMatrix [
    0.0 0.0 0.0
    1/3 2/3 0.0
    1/10 6/10 3/10
]

# ISK is guaranteed to be a an array of length R
function _weno_weights(ISk, ::Val{R}; ε = 1.0e-6, p = 2) where {R}
    C_kR = _WENO_OPTIMAL_STENCIL_WEIGHTS[:, R]
    alpha_k = map((C, IS) -> C*inv((ε+IS)^p), C_kR, ISk)
    w_k = alpha_k / sum(alpha_k)
    return w_k
end

And what I would like to do, if ISk is an AbstractArray of ForwardDiff.Dual{...}, is set the partials to zero, but I’m not sure what the recommended way to do this is. Here’s my solution so far, but I don’t think this will work with chunking.

# We want to enforce that w_k is always passive
# (to avoid information moving downwind)
function _weno_weights(
    ISk::AbstractArray{<:Dual{T}},
    order::Val{R};
    ε = 1.0e-6,
    p = 2,
) where {T,R}
    C_kR = _WENO_OPTIMAL_STENCIL_WEIGHTS[:, R]
    alpha_k = map((C, IS) -> C * inv((ε + value(IS))^p), C_kR, ISk)
    tot = sum(alpha_k)
    w_k = map(alpha_k) do a
        Dual{T}(a / tot, zero(a))
    end
    return w_k
end

Is there a recommended way to do this? Is there a convenience function/macro/example I could look at anywhere?

2

There’s ForwardDiff.value but it’s worth asking yourself whether you actually want to do that and why. In particular, it may screw up higher-order stuff if someone ever tries to differentiate through your own code