new dynamics() and output() methods in StateSpace #546
Conversation
|
Following #83 I don't see why anyone won't write def axplusbu(x, u): return A @ x + B @ uInstead what you might want to do is to provide a callable so that people can plug it in wherever they want with all bells and whistles included. In english, you return a function which is built out of the parts of the IOsystem and that is independent from the model itself. Otherwise this doesn't seem like doing anything too useful. |
|
What happens if you have a Also, what happens if a system has parameters (allowed for I/O systems)? Should I be able to set those parameter values when I call the (public) |
|
@ilayn I agree, it's not too onerous to write that out. I added these because I have been playing around with such an integrator and it does seem to clean up/improve code appearance to be able to write can you elaborate more on what you had in mind?
I see 2 options:
In
NOt sure how params are used in |
|
No matter what, I think we want the If we leave in
Another possibility would be to define However, I would leave the signatures for the signatures of I'll need to look through the FYI: for the |
…irst; more coverage in unit tests
…amics and output directly
|
OK good point, it looks like there's a pretty strong convention of leaving The latest commit changes With this, the convention now is that the call signature is |
Yes instead of doing that you can actually call the dynamics or the output methods and obtain a function so it becomes .....
def dynamics(self):
def f(x, u):
return self.a.copy() @ x + self.b.copy() @ u
return f
....And then you can ask for that function and use it elsewhere So even P is gone f will survive and become independent from the model itself hence the method becomes a constructor. But maybe this is not what you want. |
|
It seems like in the current version, the 't' argument isn't really "optional". In particular, if you don't want to specify 't' I think you would have to call the function as An alternative would be to define If you do that change, I would leave |
| dx/dt or x[t+dt] : ndarray of shape (self.nstates,) | ||
| """ | ||
|
|
||
| out = np.zeros((self.nstates, 1)) |
There was a problem hiding this comment.
I might be missing something, but vectors can have ndim=1, and the matrix products and so on will work fine. Do you need ndim=2 for the vectors here?
There was a problem hiding this comment.
you're right. I thought there might be ambiguity in how dot works on 1D vectors vs column vectors on square arrays, but it looks like it does the right thing for 1D vectors: it interprets them as column vectors. so new PR will return 1D vectors.
|
Have decided I should start over and create a new PR. In the meantime I thought I would do a quick benchmark to find out whether there are notable speed increases using a callable rather than the import control as ct
import time
sys = ct.ss([[-1.37394092, 0.08654488],
[-0.47735689, -0.43755742]], [[-0.27547913],
[ 0. ]], [[-0.03284931, -0. ]], [[0.]] )
def dynamics(self, *args):
if len(args) not in (2, 3):
raise ValueError("received"+len(args)+"args, expected 2 or 3")
if np.size(args[1]) != self.nstates:
raise ValueError("len(x) must be equal to number of states")
if len(args) == 2: # received t and x, ignore t
return self.A.dot(args[1]).reshape((-1,)) # return as row vector
else: # received t, x, and u, ignore t
#x = np.reshape(args[-2], (-1, 1)) # x must be column vector
#u = np.reshape(args[-1], (-1, 1)) # u must be column vector
if np.size(args[2]) != self.ninputs:
raise ValueError("len(u) must be equal to number of inputs")
return (self.A.dot(args[1]) + self.B.dot(args[2])).reshape((-1,)) # row
setattr(ct.StateSpace, 'dynamics', dynamics)
def dynamics_callable(sys):
return lambda t, x, u: sys.A.copy()@x + sys.B.copy()@u
dynamics = dynamics_callable(sys)
time1 = time.time()
x = sp.integrate.solve_ivp(sys.dynamics, (0,1000), (0,0), args=((1,),))
print(time.time() - time1)
time2 = time.time()
x = sp.integrate.solve_ivp(dynamics, (0,1000), (0,0), args=((1,),))
print(time.time() - time2)prints out: |
makes sense.
One possible problem arises is if To solve that ambiguity without keyword arguments, it seems like
I can do that. Though it appears that arg processing doesn't actually seem to impact speed much, see above. |
|
For |
|
I’ll add a note in the doc strings about the different anticipated use cases for |
Creates new
dynamics(x, u)andoutput(x, u)methods forStateSpacesystems. These correspond to the right hand side of the dynamcs and output equations for such systems (i.e. calculates xdot in xdot = Ax+Bu and y in y=Cx+Bu). This can be used for simulating or numerically integrating such systems in your own code, as suggested in #83.Renames equivalent private methods
_rhsand_outfunctions iniosysto same name. The latter have a slightly different call signature: the first argument is the timetrather than the statex. Since they were private methods, we can change them; it might make sense to rearrange arguments so thattis a keyword argument that appears at the end with a default value (e.g. 0) so that call signatures then are the same across system types.Follows the discussion in #434.
I considered a few alternatives for the name
dynamics:eval,evaluate,update,evolve,evolve_state,calculate, andcompute. I settled ondynamicsbecause it works for both cont-time and discrete-time systems, suggests that it is evaluating how they change in time, and doesn't suggest any mutation is happening to the system, unlikeupdate.Happy to consider alternatives.