Yes it probably shouldn't be there.

yeah i agree. Then we can simplify the API from MatrixOp(A; update_func = .., update_func! = ..), to just having one kwarg MatrixOp(A; update_func = ..), where update_func can have a (u, p, t), and a (A, u, p, t) method.

Also, I'd like to remove u dependence in update_coefficients/update_func for the reason below.

https://github.com/SciML/SciMLOperators.jl/blob/master/src/interface.jl#L42-L54

It leads to unsafe behaviour that we have no way to protect against. Removing u from update_coeffs won't lead to any depletion in use cases as u is still available at operator evaluation time.

# OOP signature
L = update_coefficients(L, p, t; kw...)
L(u, p, t)

It's very useful for defining operators used in Magnus methods though, i.e. u' = A(u)*u operations. See https://docs.sciml.ai/DiffEqDocs/stable/solvers/nonautonomous_linear_ode/#State-Dependent-Solvers . Also Jacobian-vector products and Vector-Jacobian products. This comes up in exponential Rosenbrock methods because you want to split u' = f into u' = f'u + g where g = f - f'u, and f' needs to be u-dependent. So I think it's generally a requirement to be able to have operators that are u-dependent. But I cannot see a case where a past immutable operator is used, and the update_coefficients! would only have the old operator representation to mutate it.