Re: [python-control-discuss] forced response problem for large dt

[Ryan K. <rk...@si...>]

My hackish solution to the immediate problem is to integrate 3 times at
dt/3 within the for loop:

Ndt = 3#factor to divide dt by
dt_sim = dt/Ndt

for i in range(1,N):
    e[i] = r[i] - y[i-3]#<-- one step time delay on the measurement
    #v[i] = Gc*e[i]#<-- this is only doing P control
    esum += e[i]
    v[i] = Kp*e[i] + Ki*esum + Kd*(e[i]-e[i-1]) #PID control
    if include_sat:
        v[i] = mysat(v[i])
    t0 = dt*(i-1)

    #integrate Ndt times to get to the next dt
    for q in range(Ndt):
        t1 = t0 + dt_sim
        to, yo, xo = control.forced_response(G_int, [t0,t1], [v[i],v[i]],
X0=X0)
        X0 = xo[:,-1]#<-- save for next time through for loop
        t0 = t1#<-- for next pass

    x[i] = squeeze(X0)
    y[i] = yo[-1]

This keep my students moving forward for now.


--
Ryan Krauss, Ph.D.
Associate Professor
Mechanical Engineering
Southern Illinois University Edwardsville

On Fri, Dec 4, 2015 at 3:21 PM, Ryan Krauss <rk...@si...> wrote:

> Apparently I don't know which email list is subscribed.  Sorry if this
> comes through more than once.
>
> I am teaching a classical controls class for mechanical engineering
> undergraduates.  Up till now, I have tried to gloss over state-space.  My
> students need to do an initial condition simulation for a system that
> includes actuator saturation.  I have them do this by integrating for one
> time step at a time using control.forced_response inside a for loop.  We
> are essentially doing a continuous time approximation of ZOH with the input
> held constant for each time step.  As the time step gets too large, sp.
> linalg.expm has to use a higher order pade approximation and eventually
> throws this error:
>
> ValueError                                Traceback (most recent call last)
>
>
> /Users/rkrauss/git/python-control/zumo_PID_simulation_modified_ss.py in <module>()
>
>     103     t0 = dt*(i-1)
>
>     104     t1 = dt*i
>
> --> 105
>   to, yo, xo = control.forced_response(G_int, [t0,t1], [v[i],v[i]], X0=X0)
>
>     106     X0 = xo[:,-1]#<-- save for next time through for loop
>
>     107     x[i] = squeeze(X0)
>
>
> /Users/rkrauss/git/python-control/control/timeresp.pyc in forced_response(sys,
> T, U, X0, transpose)
>
>     374                          [np.zeros((n_inputs, n_states + 2 *
> n_inputs))]])
>
>     375             print('M=' + str(M))
>
> --> 376             expM = sp.linalg.expm(M)
>
>     377             Ad = expM[:n_states, :n_states]
>
>     378             Bd1 = expM[:n_states, n_states+n_inputs:]
>
>
> /usr/local/lib/python2.7/site-packages/scipy/linalg/matfuncs.pyc in expm(A,
> q)
>
>     258     # Input checking and conversion is provided by
> sparse.linalg.expm().
>
>     259     import scipy.sparse.linalg
>
> --> 260     return scipy.sparse.linalg.expm(A)
>
>     261
>
>     262
>
>
>
> /usr/local/lib/python2.7/site-packages/scipy/sparse/linalg/matfuncs.pyc in expm(A)
>
>     580
>
>     581     """
>
> --> 582     return _expm(A, use_exact_onenorm='auto')
>
>     583
>
>     584
>
>
> /usr/local/lib/python2.7/site-packages/scipy/sparse/linalg/matfuncs.pyc in _expm(A,
> use_exact_onenorm)
>
>     635     if structure == UPPER_TRIANGULAR:
>
>     636         # Invoke Code Fragment 2.1.
>
> --> 637         X = _fragment_2_1(X, h.A, s)
>
>     638     else:
>
>     639         # X = r_13(A)^(2^s) by repeated squaring.
>
>
> /usr/local/lib/python2.7/site-packages/scipy/sparse/linalg/matfuncs.pyc in _fragment_2_1(X,
> T, s)
>
>     753     exp_diag = np.exp(scale * diag_T)
>
>     754     for k in range(n):
>
> --> 755         X[k, k] = exp_diag[k]
>
>     756
>
>     757     for i in range(s-1, -1, -1):
>
>
> ValueError: setting an array element with a sequence.
>
>
> This is probably ultimately a problem for the scipy people, but my
> students' project is due in 6 days.  Any suggestions to quickly get the
> attached simulation code to work for 60Hz simulation, i.e. dt = 1.0/60?
>
>
> Thanks,
>
>
> Ryan
>
> --
> Ryan Krauss, Ph.D.
> Associate Professor
> Mechanical Engineering
> Southern Illinois University Edwardsville
>