Let's say we have a linear differential operator $\hat A = \sum_k a_k \frac {d^k} {dx^k}$ and we want to solve the eigenvalue equation $\hat A f(x) = a f(x)$ which is an ODE that we put in sage and outputs a solution $g(x)$. We now want to substitute $f(x)$ with $g(x)$ and simplify our expression in order to extract the eigenvalues of $\hat A$. The problem here is that I am not able to substitute the derivatives with my solution neither via eqn.substitute_expression(f(x) == g(x)) nor via eqn.substitute_function(f(x),g(x)) because D[0]f(x), ... D[0,0,...,0]f(x) remain unchanged.