Scalar / Array mismatches in iosys not compatible with NumPy 1.24 #853

Activity

1. bnavigator commented on Jan 26, 2023

   

   bnavigator

   on Jan 26, 2023

   ContributorAuthor

   More actions

   When executing the example:

   In [4]: eqpt = ct.find_eqpt(io_predprey, X0, 0)
      ...: xeq = eqpt[0]                         # choose the nonzero equilibrium point
      ...: lin_predprey = ct.linearize(io_predprey, xeq, 0)
   ---------------------------------------------------------------------------
   ValueError                                Traceback (most recent call last)
   Cell In[4], line 3
         1 eqpt = ct.find_eqpt(io_predprey, X0, 0)
         2 xeq = eqpt[0]                         # choose the nonzero equilibrium point
   ----> 3 lin_predprey = ct.linearize(io_predprey, xeq, 0)
   
   File /usr/lib/python3.10/site-packages/control/iosys.py:2252, in linearize(sys, xeq, ueq, t, params, **kw)
      2250 if not isinstance(sys, InputOutputSystem):
      2251     raise TypeError("Can only linearize InputOutputSystem types")
   -> 2252 return sys.linearize(xeq, ueq, t=t, params=params, **kw)
   
   File /usr/lib/python3.10/site-packages/control/iosys.py:580, in InputOutputSystem.linearize(self, x0, u0, t, params, eps, name, copy_names, **kwargs)
       578     du = np.zeros((ninputs,))
       579     du[i] = eps
   --> 580     B[:, i] = (self._rhs(t, x0, u0 + du) - F0) / eps
       581     D[:, i] = (self._out(t, x0, u0 + du) - H0) / eps
       583 # Create the state space system
   
   File /usr/lib/python3.10/site-packages/control/iosys.py:866, in NonlinearIOSystem._rhs(self, t, x, u)
       863 def _rhs(self, t, x, u):
       864     xdot = self.updfcn(t, x, u, self._current_params) \
       865         if self.updfcn is not None else []
   --> 866     return np.array(xdot).reshape((-1,))
   
   ValueError: setting an array element with a sequence. The requested array has an inhomogeneous shape after 1 dimensions. The detected shape was (2,) + inhomogeneous part.
2. bnavigator

   changed the title [-]Avoid ragged shapes in linearize()[/-] [+]Avoid ragged shapes in linearize() / iosys[/+] on Jan 26, 2023
3. bnavigator commented on Jan 26, 2023

   

   bnavigator

   on Jan 26, 2023

   ContributorAuthor

   More actions

   This would fix the error:

   def predprey_rhs(t, x, u, params):
       # Parameter setup
       a = params.get('a', 3.2)
       b = params.get('b', 0.6)
       c = params.get('c', 50.)
       d = params.get('d', 0.56)
       k = params.get('k', 125)
       r = params.get('r', 1.6)
   
       # Map the states into local variable names
       H = x[0]
       L = x[1]
   
       # Compute the control action (only allow addition of food)
       u_0 = u[0] if u[0] > 0 else 0
   
       # Compute the discrete updates
       dH = (r + u_0) * H * (1 - H/k) - (a * H * L)/(c + H)
       dL = b * (a * H *  L)/(c + H) - d * L
   
       return np.array([dH, dL])

   Documentation says, that the update function's u must be an array. As well as the ueq for linearize:

   python-control/control/iosys.py

   Lines 2198 to 2214 in 55730b6

   	def linearize(sys, xeq, ueq=None, t=0, params=None, **kw):
   	"""Linearize an input/output system at a given state and input.
   	This function computes the linearization of an input/output system at a
   	given state and input value and returns a :class:`~control.StateSpace`
   	object. The evaluation point need not be an equilibrium point.
   	Parameters
   	----------
   	sys : InputOutputSystem
   	The system to be linearized
   	xeq : array
   	The state at which the linearization will be evaluated (does not need
   	to be an equilibrium state).
   	ueq : array
   	The input at which the linearization will be evaluated (does not need
   	to correspond to an equlibrium state).

   Strictly speaking, the zero us from the example would not be allowed:

   eqpt = ct.find_eqpt(io_predprey, X0, 0)
   xeq = eqpt[0]                         # choose the nonzero equilibrium point
   lin_predprey = ct.linearize(io_predprey, xeq, 0)

   But we have this:

   python-control/control/iosys.py

   Lines 547 to 551 in 55730b6

   	# Convert x0, u0 to arrays, if needed
   	if np.isscalar(x0):
   	x0 = np.ones((nstates,)) * x0
   	if np.isscalar(u0):
   	u0 = np.ones((ninputs,)) * u0
4. bnavigator

   changed the title [-]Avoid ragged shapes in linearize() / iosys[/-] [+]Scalar / Array mismatches in iosys not compatible with NumPy 1.24[/+] on Jan 26, 2023
5. bnavigator commented on Jan 26, 2023

   

   bnavigator

   on Jan 26, 2023

   Last edited by bnavigator

   ContributorAuthor

   More actions

   The test suite has a different predprey model:

   python-control/control/tests/iosys_test.py

   Lines 1558 to 1571 in c3488cd

   	def predprey(t, x, u, params={}):
   	"""Predator prey dynamics"""
   	r = params.get('r', 2)
   	d = params.get('d', 0.7)
   	b = params.get('b', 0.3)
   	k = params.get('k', 10)
   	a = params.get('a', 8)
   	c = params.get('c', 4)
   	# Dynamics for the system
   	dx0 = r * x[0] * (1 - x[0]/k) - a * x[1] * x[0]/(c + x[0])
   	dx1 = b * a * x[1] * x[0] / (c + x[0]) - d * x[1]
   	return np.array([dx0, dx1])
6. murrayrm

   closed this as completedin #854on Jan 27, 2023