scipy ode update set_f_params inside function set as set_solout

When integrating an ode with scipy, ode accepts a function with more arguments than t and y. For example:

def fun(t, y, param1, param2):

and the value of these arguments can be set using set_f_params method.

However, when using also set_solout method and trying to update the params with set_f_params inside this function, the integration remains the same as if the params were not being modified.

How would you modify the the params using sol_out? I would like to benefit from dopri5 dense output, but I need the non-homogeneous terms to be updated at every time step.

A minimal example is shown below.

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import ode

def fun(t, x, param):
    return x - param

def f_param(t):
    return t

ode1 = ode(fun).set_integrator('dopri5').set_initial_value([10.0])
ode1.set_f_params(f_param(0))
results1 = ([], [])

ode2 = ode(fun).set_integrator('dopri5').set_initial_value([10.0])
ode2.set_f_params(f_param(0))
results2 = ([], [])

def callback1(t, x):
    results1[0].append(t)
    results1[1].append(x.copy())

def callback2(t, x):
    results2[0].append(t)
    results2[1].append(x.copy())
    ode2.set_f_params(f_param(t))

ode1.set_solout(callback1)
ode2.set_solout(callback2)

ode1.integrate(3)
ode2.integrate(3)

plt.plot(results1[0], results1[1], 'o-', alpha=0.7, label='ode1')
plt.plot(results2[0], results2[1], '.--', label='ode2')
plt.legend()

and the results are shown here: