I’m moving my first steps in programming with Julia and I have written a module whose structure is below:
module my_solver
function solve(f,u0)
# use f(u) and u0 to compute u_final
return u_final
end
end
in which f is specified so far as a non in-place function. I am now in the process to update my module so that it can deal with both in place and not in place functions. At fist I would think of modifying the module as follows
module my_solver
function solve_not_in_place(f, u0)
# use f(u) and u0 to compute u_final
return u_final
end
function solve_in_place(f!, u0)
z = similar(u0)
# use f!(z, u) and u0 to compute u_final
return u_final
end
end
Questions I have:
I would like to provide my users with a single method, say solve, that understands the appropriate input and calls the appropriate function. What’s an appropriate Julian pattern for that?
Might it be that I can define two functions solve(f,u0) and let Julia sort out whether the function is in place or not?
I have assumed above that the best practice here is to split the two functions, as opposed to merge them. Is that advisable from your viewpoint (I understand that the answer may be dependent on what gets on each function)?
Of course any other comments or suggestions are greatly appreciated
Hard to say from the pseudocode because we don’t see what you do with the f and f! calls or how u_final gets defined. So far it’s not obvious exactly where the problem would be that justifies splitting into 2 solve_X functions.
Speaking more generally, in-place methods often return the key mutated object (so f!(key, args...) returns key) so they can work in the same return value assignments (result = f_(key, args...)) that work for non-in-place methods. If there isn’t one key mutated object to return or it doesn’t make sense to return all of the mutated objects, then it just shouldn’t be an input for the higher-order functions that take non-in-place methods.
It’s tempting to quote Feyerabend, “anything goes”. Although it might be possible, it’s a bit non-julian to try to figure out if the passed function is in-place. A more julian way is to wrap the function:
struct InPlace{FT} <: Function
f::FT
end
(fs::InPlace)(x...) = fs.f(x...)
solve(f::InPlace, x0) = inplacesolver(f, x0)
solve(f, x0) = outofplacesolver(f, x0)
Such a wrapper may be made mandatory, and you may wrap f together with x0. It all depends on typical usage.
However, a better place to start is to adopt the CommonSolve.jl framework. So that your solver looks like other solvers.
I agree that in-place functions should always be accomodated. A similar thing is done with LinearAlgebra.mul!, which computes a*A*B + b*C for matrices A,B,C, scalars a,b, in-place in C.
It’s a matter of taste if you implement the in-place solver as solve!(f, u0), or solve(InPlace(f), u0), or solve(InPlace(f, u0)) as I suggested above.
Sorry I should have been clearer. Your comments already more or less communicate what function calls are being used. Pseudocode communicates more than that, like what happens to the variables immediately around the calls. For example, that 2nd bullet point could be u_final = f(u) or it could be u_final = u; f!(u_final).
That said, a non-in-place f would require u_final = f(u), and an in-place f! can be used in the same expression u_final = f!(u) if it returns its mutated input.
Good to see you on here! Hope to see you doing lots of things in Julia I’d recommend looking at the code in SciMLBase.jl, specifically for the ODEProblem (the other structures behave similarly). There a problem object is created similar to sgaure’s reply. There is also some autodetection of whether the function is inplace or not.
For many of the things you work on, there is probably a fair bit of inspiration that can be taken from the SciML ecosystem.
I’d recommend looking at the code in SciMLBase.jl, specifically for the