RFC: storing intermediate values
Sometimes I solve nonlinear systems of the form
f(g(x)) ≈ 0 # could be a scalar or a vector
where g(x) is a large, costly, and ultimately interesting object, so one would call it on the root, but that is a bit wasteful. It would be nice to have an interface that returns it. I am thinking of something like
(root = x, value = g(x)) = find_zero(f, arg_or_bracket, method; by = identity)
where by would be g, and a NamedTuple is returned. But this of course would be breaking, so a new name would perhaps be better.
Why can't you just use a closure or something and store it yourself?
local val
find_zero(x -> (val = g(x); f(val)))
I don't know how to get the type inferred correctly. Eg
using Roots
function find_zero_ext(f, g, bracket)
local val
x = find_zero(x -> (val = g(x); f(val)), bracket)
(x, val)
end
results in
julia> @code_warntype find_zero_ext(x -> x^3, cbrt, (-1, 1))
Body::Tuple{Float64,Any}
1 ─ %1 = %new(Core.Box)::Core.Box
│ %2 = %new(getfield(Main, Symbol("##9#10")){getfield(Main, Symbol("##15#16")),typeof(cbrt)}, f, g, %1)::getfield(Main, Symbol("##9#10")){getfield(Main, Symbol("##15#16")),typeof(cbrt)}
│ %3 = Roots.find_zero::typeof(find_zero)
│ %4 = invoke Roots.:(#find_zero#8)($(QuoteNode(Roots.NullTracks()))::Roots.NullTracks, false::Bool, $(QuoteNode(Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}()))::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, %3::Function, %2::Function, _4::Tuple{Int64,Int64}, $(QuoteNode(Bisection()))::Bisection)::Float64
│ %5 = (Core.isdefined)(%1, :contents)::Bool
└── goto #3 if not %5
2 ─ goto #4
3 ─ $(Expr(:throw_undef_if_not, :val, false))
4 ┄ %9 = (Core.getfield)(%1, :contents)::Any
│ %10 = (Core.tuple)(%4, %9)::Tuple{Float64,Any}
└── return %10
That's the old "boxes in closures"-issue. Can be worked around with something like:
julia> function find_zero_ext(f, g, bracket)
val = Ref{Float64}()
x = find_zero(x -> (val[] = g(x); f(val[])), bracket)
(x, val[])
end
Assume I don't know the result type of g(x) (or that it's a nested NamedTuple of epic horribleness).