Hello, I would like to solve the ODE dy/dt = -2y + data(t), between t=0..4, for y(t=0)=1. I
wrote the following code: import numpy as np from scipy.integrate import odeint from scipy.interpolate import interp1d t = np.linspace(0, 3, 4) data = [1, 2, 3, 4] linear_interpolation = interp1d(t, data) def func(y, t0): print 't0', t0 return -2*y + linear_interpolation(t0) soln = odeint(func, 1, t) When I run this code, I get several errors: ValueError: A value in x_new is above the interpolation range. odepack.error: Error occurred while calling the Python function named func My interpolation range is between 0.0 and 3.0. Printing the value of t0 in func, I realized that t0 is actually sometimes above my interpolation range: 3.07634612585, 3.0203768998, 3.00638459329, ... I have a few questions: - how does integrate.ode makes t0 vary? Why does it make t0 exceed the infimum (3.0) of my interpolation range? - in spite of these errors, integrate.ode returns an array which seems to contain correct value. So, should I just catch and ignore these
errors? - if I shouldn't ignore these errors, what is the best way to avoid them? 2 suggestions for the last question: - in interp1d, I could set bounds_error=False and fill_value=data[-1] since the t0 outside of my interpolation range seem to be closed to t[-1]: linear_interpolation = interp1d(t, data, bounds_error=False, fill_value=data[-1]) But first I would like to be sure that with any other func and any other data the t0 will always remain closed to t[-1]. For example, if integrate.ode chooses a t0 below my interpolation range, the fill_value would be still data[-1], which would not be correct. Maybe to know how integrate.ode makes t0 vary would help me to be sure of that (see my first question). - in func, I could enclose the linear_interpolation call in a try/except block, and, when I catch a ValueError, I recall linear_interpolation but with t0 truncated: def func(y, t0): try: interpolated_value = linear_interpolation(t0) except ValueError: interpolated_value = linear_interpolation(int(t0)) # truncate t0 return -2*y + interpolated_value At least this solution permits linear_interpolation to still raise an exception if integrate.ode makes a t0 above 4.0 or below -1.0. I can then be alerted of incoherent behavior. But it is not really readable and the truncation seems to me a little arbitrary by now. Maybe I'm just overthinking about these errors. Please let me know. Thanks in advance. Cheers, Camille _______________________________________________ SciPy-User mailing list [hidden email] http://mail.scipy.org/mailman/listinfo/scipy-user |
camille chambon <camillechambon <at> yahoo.fr> writes: > > > Hello, > > I would like to solve the ODE dy/dt = -2y + data(t), between t=0..3, for y(t=0)=1. > > I > wrote the following code: > > import numpy as npfrom scipy.integrate import odeintfrom scipy.interpolate import interp1dt = np.linspace(0, 3, 4)data = [1, 2, 3, 4]linear_interpolation = interp1d(t, data)def func(y, t0): > print 't0', t0 return -2*y + linear_interpolation(t0)soln = odeint(func, 1, t) > > When I run this code, I get several errors: > > > ValueError: A value in x_new is above the interpolation range.odepack.error: Error occurred while calling the Python function named func > > My interpolation range is between 0.0 and 3.0. > > Printing the value of t0 in func, I realized that t0 is actually sometimes above my interpolation range: 3.07634612585, 3.0203768998, 3.00638459329, ... > > I have a few questions: > > - how does integrate.ode makes t0 vary? Why does it make t0 exceed the infimum (3.0) of my interpolation range? > > - in spite of these errors, integrate.ode returns an array which seems to contain correct value. So, should I just catch and ignore these > errors? > > > - if I shouldn't ignore these errors, what is the best way to avoid them? > > > 2 suggestions for the last question: > > - in interp1d, I could set bounds_error=False and fill_value=data[-1] since the t0 outside of my interpolation range seem to be closed to t[-1]: > linear_interpolation = interp1d(t, data, bounds_error=False, fill_value=data[-1]) > But first I would like to be sure that with any other func and any other data the t0 will always remain closed to t[-1]. For example, if integrate.ode chooses a t0 below my interpolation range, the fill_value would be still data[-1], which would not be correct. Maybe to know how integrate.ode makes t0 vary would help me to be sure of that (see my first question). > > - in func, I could enclose the linear_interpolation call in a try/except block, and, when I catch a ValueError, I recall linear_interpolation but with t0 truncated: > > def func(y, t0): > try: > interpolated_value = linear_interpolation(t0) > except ValueError: > > interpolated_value = linear_interpolation(int(t0)) # truncate t0 > > > return -2*y + interpolated_value > > > > > > At least this solution permits linear_interpolation to still raise an alerted of incoherent behavior. But it is not really readable and the truncation seems to me a little arbitrary by now. > > > Maybe I'm just overthinking about these errors. Please let me know. > > Thanks in advance. > > Cheers, > > Camille > > > > _______________________________________________ > SciPy-User mailing list > SciPy-User <at> scipy.org > http://mail.scipy.org/mailman/listinfo/scipy-user > Sorry, I made a typo in the range of t: I want to solve the ODE dy/dt = -2y + data(t) between t=0..3, and not between t=0..4. I corrected my original message. Cheers, Camille _______________________________________________ SciPy-User mailing list [hidden email] http://mail.scipy.org/mailman/listinfo/scipy-user |
Salut Camille:
Just don't ask for the solution in the last point, since in order to calculate a value of the solution at that point it has to go a bit beyond (y should be defined from both sides of ti for the derivative to exist in ti (this is a necessary condition)), but interpolation does not really define values beyond the right limit. Here I've modified your example: http://nbviewer.ipython.org/github/guziy/PyNotebooks/blob/master/ode_demo.ipynb Cheers 2014-07-30 4:05 GMT-04:00 camille chambon <[hidden email]>: > > camille chambon <camillechambon <at> yahoo.fr> writes: > >> >> >> Hello, >> >> I would like to solve the ODE dy/dt = -2y + data(t), between t=0..3, for > y(t=0)=1. >> >> I >> wrote the following code: >> >> import numpy as npfrom scipy.integrate import odeintfrom scipy.interpolate > import interp1dt = np.linspace(0, 3, 4)data = [1, 2, 3, > 4]linear_interpolation = interp1d(t, data)def func(y, t0): >> print 't0', t0 return -2*y + linear_interpolation(t0)soln = > odeint(func, 1, t) >> >> When I run this code, I get several errors: >> >> >> ValueError: A value in x_new is above the interpolation > range.odepack.error: Error occurred while calling the Python function named func >> >> My interpolation range is between 0.0 and 3.0. >> >> Printing the value of t0 in func, I realized that t0 is actually sometimes > above my interpolation range: 3.07634612585, 3.0203768998, 3.00638459329, ... >> >> I have a few questions: >> >> - how does integrate.ode makes t0 vary? Why does it make t0 exceed the > infimum (3.0) of my interpolation range? >> >> - in spite of these errors, integrate.ode returns an array which seems to > contain correct value. So, should I just catch and ignore these >> errors? >> >> >> - if I shouldn't ignore these errors, what is the best way to avoid them? >> >> >> 2 suggestions for the last question: >> >> - in interp1d, I could set bounds_error=False and fill_value=data[-1] > since the t0 outside of my interpolation range seem to be closed to t[-1]: >> linear_interpolation = interp1d(t, data, bounds_error=False, > fill_value=data[-1]) >> But first I would like to be sure that with any other func and any other > data the t0 will always remain closed to t[-1]. For example, if > integrate.ode chooses a t0 below my interpolation range, the fill_value > would be still data[-1], which would not be correct. Maybe to know how > integrate.ode makes t0 vary would help me to be sure of that (see my first > question). >> >> - in func, I could enclose the linear_interpolation call in a try/except > block, and, when I catch a ValueError, I recall linear_interpolation but > with t0 truncated: >> >> def func(y, t0): >> try: >> interpolated_value = linear_interpolation(t0) >> except ValueError: >> >> interpolated_value = linear_interpolation(int(t0)) # truncate t0 >> >> >> return -2*y + interpolated_value >> >> >> >> >> >> At least this solution permits linear_interpolation to still raise an > exception if integrate.ode makes a t0 above 4.0 or below -1.0. I can then be > alerted of incoherent behavior. But it is not really readable and the > truncation seems to me a little arbitrary by now. >> >> >> Maybe I'm just overthinking about these errors. Please let me know. >> >> Thanks in advance. >> >> Cheers, >> >> Camille >> >> >> >> _______________________________________________ >> SciPy-User mailing list >> SciPy-User <at> scipy.org >> http://mail.scipy.org/mailman/listinfo/scipy-user >> > > Sorry, I made a typo in the range of t: I want to solve the ODE dy/dt = -2y > + data(t) between t=0..3, and not between t=0..4. I corrected my original > message. > Cheers, > Camille > > > _______________________________________________ > SciPy-User mailing list > [hidden email] > http://mail.scipy.org/mailman/listinfo/scipy-user -- Sasha _______________________________________________ SciPy-User mailing list [hidden email] http://mail.scipy.org/mailman/listinfo/scipy-user |
In reply to this post by CamilleC
On Wed, Jul 30, 2014 at 4:05 AM, camille chambon <[hidden email]> wrote: > I would like to solve the ODE dy/dt = -2y + data(t), between t=0..3, for Camille, I added an answer to your question on StackOverflow: http://stackoverflow.com/questions/25031966/integrate-ode-sets-t0-values-outside-of-my-data-range/ Summary: * It is normal for the ODE solver to evaluate your function at points beyond the last requested time value. * One way to avoid the interpolation problem is to extend the interpolation data linearly, using the last two data point. Warren
_______________________________________________ SciPy-User mailing list [hidden email] http://mail.scipy.org/mailman/listinfo/scipy-user |
In reply to this post by Oleksandr Huziy
Oleksandr Huziy <guziy.sasha <at> gmail.com> writes:
> > Salut Camille: > > Just don't ask for the solution in the last point, since in order to > calculate a value of the solution at that point it has to go a bit > beyond (y should be defined from both sides of ti for the derivative > to exist in ti (this is a necessary condition)), but interpolation > does not really define values beyond the right limit. > > Here I've modified your example: > > > Cheers > Hello Oleksandr, Thanks for your answer. But I do need the solution in the last point. Furthermore, in your modifications you redefine data(t) as: [1, 2, 3, 4]=data([0.,1.33333333,2.66666667,4.], which is not the same as: [1, 2, 3, 4]=data([0.,1.,2.,3.]. I need data(t) to remain as: [1, 2, 3, 4]=data([0.,1.,2.,3.]. Cheers, Camille _______________________________________________ SciPy-User mailing list [hidden email] http://mail.scipy.org/mailman/listinfo/scipy-user |
In reply to this post by Warren Weckesser-2
Warren Weckesser <warren.weckesser <at> gmail.com> writes:
> > Camille, > > I added an answer to your question on StackOverflow: http://stackoverflow.com/questions/25031966/integrate-ode-sets-t0-values-outside-of-my-data-range/ > > Summary: > > * It is normal for the ODE solver to evaluate your function at points beyond the last requested time value. > > * One way to avoid the interpolation problem is to extend the interpolation data linearly, using the last two data point. > > Warren > > Hello Warren, Thanks for your answer. Actually, data(t) is a sinusoidal function. Thus I can not extend the interpolation data linearly. I edited my question on StackOverflow accordingly. Cheers, Camille _______________________________________________ SciPy-User mailing list [hidden email] http://mail.scipy.org/mailman/listinfo/scipy-user |
On 7 August 2014 11:05, Camille <[hidden email]> wrote: Thanks for your answer. Actually, data(t) is a sinusoidal function. Thus I Note that you are evaluating it very close to the boundary, so the interpolation effects will not be so bad. You can check the sensibility comparing the results with a purposefully "wrong" interpolation, like the same linear interpolation but with the opposite slope; but I bet the differences are going to be slim. If you know your data is sinusoidal, you can use that to make an even better estimation of the next value. Essentially, you need to provide the ODE solver with a way to estimate the derivatives of your function at any point. /David _______________________________________________ SciPy-User mailing list [hidden email] http://mail.scipy.org/mailman/listinfo/scipy-user |
On 8/7/14, Daπid <[hidden email]> wrote:
> On 7 August 2014 11:05, Camille <[hidden email]> wrote: > >> Thanks for your answer. Actually, data(t) is a sinusoidal function. Thus >> I >> can not extend the interpolation data linearly. I edited my question on >> StackOverflow accordingly. > > > Note that you are evaluating it very close to the boundary, so the > interpolation effects will not be so bad. You can check the sensibility > comparing the results with a purposefully "wrong" interpolation, like the > same linear interpolation but with the opposite slope; but I bet the > differences are going to be slim. > > If you know your data is sinusoidal, you can use that to make an even > better estimation of the next value. Essentially, you need to provide the > ODE solver with a way to estimate the derivatives of your function at any > point. > > > /David > I agree with David. Use whatever extrapolation method is appropriate for your function. The main point is to expect your function to be evaluated a little bit beyond the final time requested in odeint. I updated my answer on StackOverflow with a suggestion to use either the "dopri5" or "dop853" solver of the scipy.integrate.ode class. It appears that these solvers do not evaluate your function at times beyond the requested time. Check out the SO answer for the sample code: http://stackoverflow.com/questions/25031966/integrate-ode-sets-t0-values-outside-of-my-data-range Warren _______________________________________________ SciPy-User mailing list [hidden email] http://mail.scipy.org/mailman/listinfo/scipy-user |
Thanks Warren and David.
The function to be integrated might not be sinusoidal all the time and it's not possible to find an extrapolation method which is valid all the time. Thus I will use the "dopri5" or "dop853" solver of the scipy.integrate.ode class. Thanks again. Camille |
Free forum by Nabble | Edit this page |