[SciPy-User] bug in signal.lsim2

Warren Weckesser warren.weckesser@enthought....
Thu Jan 28 20:05:28 CST 2010


Or, you could modify your script lsim2_problem.py to eliminate the 
initial interval during which the solution is exactly zero, by changing 
this line:

    u2[50:100] = amp

to this:

    u2[:50] = amp

This just translates the pulse so that it starts at t=0. Since the 
default initial condition used by lsim2 is zero, it will give the same 
solution, just translated in time.

Warren

Warren Weckesser wrote:
> Ryan,
>
> The problem is that the ODE solver used by lsim2 is too good. :)
>
> It uses scipy.integrate.odeint, which in turn uses the Fortran library 
> LSODA.  Like any good solver, LSODA is an adaptive solver--it adjusts 
> its step size to be as large as possible while keeping estimates of 
> the error bounded.  For the problem you are solving, with initial 
> condition 0, the exact solution is initially exactly 0.  This is such 
> a nice smooth solution that the solver's step size quickly grows--so 
> big, in fact, that it skips right over your pulse and never sees it.
>
> So how does it create all those intermediate points at the requested 
> time values?  It uses interpolation between the steps that it computed 
> to create the solution values at the times that you requested.  So 
> using a finer grid of time values won't help.  (If lsim2 gave you a 
> hook into the parameters passed to odeint, you could set odeint's 
> 'hmax' to a value smaller than your pulse width, which would force the 
> solver to see the pulse.  But there is no way to set that parameter 
> from lsim2.)
>
> The basic problem is you are passing in a discontinuous function to a 
> solver that expects a smooth function.  A better way to solve this 
> problem is to explicitly account for the discontinuity. One 
> possibility is the attached script.
>
> This is an excellent "learning opportunity" for your students on the 
> hazards of numerical computing!
>
> Warren
>
>
> Ryan Krauss wrote:
>> I believe I have discovered a bug in signal.lsim2.  I believe the
>> short attached script illustrates the problem.  I was trying to
>> predict the response of a transfer function with a pure integrator:
>>
>>              g
>> G = -------------
>>          s(s+p)
>>
>> to a finite width pulse.  lsim2 seems to handle the step response just
>> fine, but says that the pulse response is exactly 0.0 for the entire
>> time of the simulation.  Obviously, this isn't the right answer.
>>
>> I am running scipy 0.7.0 and numpy 1.2.1 on Ubuntu 9.04, but I also
>> have the same problem on Windows running 0.7.1 and 1.4.0.
>>
>> Thanks,
>>
>> Ryan
>>   
>> ------------------------------------------------------------------------
>>
>> _______________________________________________
>> SciPy-User mailing list
>> SciPy-User@scipy.org
>> http://mail.scipy.org/mailman/listinfo/scipy-user
>
> ------------------------------------------------------------------------
>
> ------------------------------------------------------------------------
>
> _______________________________________________
> SciPy-User mailing list
> SciPy-User@scipy.org
> http://mail.scipy.org/mailman/listinfo/scipy-user



More information about the SciPy-User mailing list