replace ndarray.dot() with @

control/canonical.py

@@ -8,7 +8,7 @@
 
 import numpy as np
 
-from numpy import zeros, zeros_like, shape, poly, iscomplex, vstack, hstack, dot, \
+from numpy import zeros, zeros_like, shape, poly, iscomplex, vstack, hstack, \
     transpose, empty, finfo, float64
 from numpy.linalg import solve, matrix_rank, eig
 
@@ -149,7 +149,7 @@ def observable_form(xsys):
         raise ValueError("Transformation matrix singular to working precision.")
 
     # Finally, compute the output matrix
-    zsys.B = Tzx.dot(xsys.B)
+    zsys.B = Tzx @ xsys.B
 
     return zsys, Tzx
 
@@ -189,13 +189,13 @@ def rsolve(M, y):
 
     # Update the system matrices
     if not inverse:
-        zsys.A = rsolve(T, dot(T, zsys.A)) / timescale
-        zsys.B = dot(T, zsys.B) / timescale
+        zsys.A = rsolve(T, T @ zsys.A) / timescale
+        zsys.B = T @ zsys.B / timescale
         zsys.C = rsolve(T, zsys.C)
     else:
-        zsys.A = solve(T, zsys.A).dot(T) / timescale
+        zsys.A = solve(T, zsys.A) @ T / timescale
         zsys.B = solve(T, zsys.B) / timescale
-        zsys.C = zsys.C.dot(T)
+        zsys.C = zsys.C @ T
 
     return zsys
 
@@ -405,8 +405,8 @@ def bdschur(a, condmax=None, sort=None):
         permidx = np.hstack([blkidxs[i] for i in sortidx])
         rperm = np.eye(amodal.shape[0])[permidx]
 
-        tmodal = tmodal.dot(rperm)
-        amodal = rperm.dot(amodal).dot(rperm.T)
+        tmodal = tmodal @ rperm
+        amodal = rperm @ amodal @ rperm.T
         blksizes = blksizes[sortidx]
 
     elif sort is None:

control/frdata.py

@@ -47,7 +47,7 @@
 from warnings import warn
 import numpy as np
 from numpy import angle, array, empty, ones, \
-    real, imag, absolute, eye, linalg, where, dot, sort
+    real, imag, absolute, eye, linalg, where, sort
 from scipy.interpolate import splprep, splev
 from .lti import LTI, _process_frequency_response
 from . import config
@@ -301,7 +301,7 @@ def __mul__(self, other):
         fresp = empty((outputs, inputs, len(self.omega)),
                       dtype=self.fresp.dtype)
         for i in range(len(self.omega)):
-            fresp[:, :, i] = dot(self.fresp[:, :, i], other.fresp[:, :, i])
+            fresp[:, :, i] = self.fresp[:, :, i] @ other.fresp[:, :, i]
         return FRD(fresp, self.omega,
                    smooth=(self.ifunc is not None) and
                           (other.ifunc is not None))
@@ -329,7 +329,7 @@ def __rmul__(self, other):
         fresp = empty((outputs, inputs, len(self.omega)),
                       dtype=self.fresp.dtype)
         for i in range(len(self.omega)):
-            fresp[:, :, i] = dot(other.fresp[:, :, i], self.fresp[:, :, i])
+            fresp[:, :, i] = other.fresp[:, :, i] @ self.fresp[:, :, i]
         return FRD(fresp, self.omega,
                    smooth=(self.ifunc is not None) and
                           (other.ifunc is not None))

control/mateqn.py

@@ -35,7 +35,7 @@
 # OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
 # SUCH DAMAGE.
 
-from numpy import shape, size, asarray, copy, zeros, eye, dot, \
+from numpy import shape, size, asarray, copy, zeros, eye, \
     finfo, inexact, atleast_2d
 from scipy.linalg import eigvals, solve_discrete_are, solve
 from .exception import ControlSlycot, ControlArgument
@@ -624,9 +624,9 @@ def care(A, B, Q, R=None, S=None, E=None, stabilizing=True):
 
         # Calculate the gain matrix G
         if size(R_b) == 1:
-            G = dot(dot(1/(R_ba), asarray(B_ba).T), X)
+            G = 1/(R_ba) * asarray(B_ba).T @ X
         else:
-            G = dot(solve(R_ba, asarray(B_ba).T), X)
+            G = solve(R_ba, asarray(B_ba).T) @ X
 
         # Return the solution X, the closed-loop eigenvalues L and
         # the gain matrix G
@@ -732,9 +732,9 @@ def care(A, B, Q, R=None, S=None, E=None, stabilizing=True):
 
         # Calculate the gain matrix G
         if size(R_b) == 1:
-            G = dot(1/(R_b), dot(asarray(B_b).T, dot(X, E_b)) + asarray(S_b).T)
+            G = 1/(R_b) * (asarray(B_b).T @ X @ E_b + asarray(S_b).T)
         else:
-            G = solve(R_b, dot(asarray(B_b).T, dot(X, E_b)) + asarray(S_b).T)
+            G = solve(R_b, asarray(B_b).T @ X @ E_b + asarray(S_b).T)
 
         # Return the solution X, the closed-loop eigenvalues L and
         # the gain matrix G
@@ -794,8 +794,8 @@ def dare(A, B, Q, R, S=None, E=None, stabilizing=True):
         Rmat = _ssmatrix(R)
         Qmat = _ssmatrix(Q)
         X = solve_discrete_are(A, B, Qmat, Rmat)
-        G = solve(B.T.dot(X).dot(B) + Rmat, B.T.dot(X).dot(A))
-        L = eigvals(A - B.dot(G))
+        G = solve(B.T @ X @ B + Rmat, B.T @ X @ A)
+        L = eigvals(A - B @ G)
         return _ssmatrix(X), L, _ssmatrix(G)
 
 
@@ -926,11 +926,11 @@ def dare_old(A, B, Q, R, S=None, E=None, stabilizing=True):
 
         # Calculate the gain matrix G
         if size(R_b) == 1:
-            G = dot(1/(dot(asarray(B_ba).T, dot(X, B_ba)) + R_ba),
-                    dot(asarray(B_ba).T, dot(X, A_ba)))
+            G = (1/(asarray(B_ba).T @ X @ B_ba + R_ba) *
+                 asarray(B_ba).T @ X @ A_ba)
         else:
-            G = solve(dot(asarray(B_ba).T, dot(X, B_ba)) + R_ba,
-                      dot(asarray(B_ba).T, dot(X, A_ba)))
+            G = solve(asarray(B_ba).T @ X @ B_ba + R_ba,
+                      asarray(B_ba).T @ X @ A_ba)
 
         # Return the solution X, the closed-loop eigenvalues L and
         # the gain matrix G
@@ -1036,11 +1036,11 @@ def dare_old(A, B, Q, R, S=None, E=None, stabilizing=True):
 
         # Calculate the gain matrix G
         if size(R_b) == 1:
-            G = dot(1/(dot(asarray(B_b).T, dot(X, B_b)) + R_b),
-                    dot(asarray(B_b).T, dot(X, A_b)) + asarray(S_b).T)
+            G = (1/(asarray(B_b).T @ X @ B_b + R_b) *
+                 (asarray(B_b).T @ X @ A_b + asarray(S_b).T))
         else:
-            G = solve(dot(asarray(B_b).T, dot(X, B_b)) + R_b,
-                      dot(asarray(B_b).T, dot(X, A_b)) + asarray(S_b).T)
+            G = solve(asarray(B_b).T @ X @ B_b + R_b,
+                      asarray(B_b).T @ X @ A_b + asarray(S_b).T)
 
         # Return the solution X, the closed-loop eigenvalues L and
         # the gain matrix G

control/statesp.py

@@ -49,8 +49,8 @@
 
 import math
 import numpy as np
-from numpy import any, array, asarray, concatenate, cos, delete, \
-    dot, empty, exp, eye, isinf, ones, pad, sin, zeros, squeeze, pi
+from numpy import any, asarray, concatenate, cos, delete, \
+    empty, exp, eye, isinf, ones, pad, sin, zeros, squeeze
 from numpy.random import rand, randn
 from numpy.linalg import solve, eigvals, matrix_rank
 from numpy.linalg.linalg import LinAlgError
@@ -1126,23 +1126,23 @@ def lft(self, other, nu=-1, ny=-1):
         H22 = TH[ny:, self.nstates + other.nstates + self.ninputs - nu:]
 
         Ares = np.block([
-            [A + B2.dot(T21), B2.dot(T22)],
-            [Bbar1.dot(T11), Abar + Bbar1.dot(T12)]
+            [A + B2 @ T21, B2 @ T22],
+            [Bbar1 @ T11, Abar + Bbar1 @ T12]
         ])
 
         Bres = np.block([
-            [B1 + B2.dot(H21), B2.dot(H22)],
-            [Bbar1.dot(H11), Bbar2 + Bbar1.dot(H12)]
+            [B1 + B2 @ H21, B2 @ H22],
+            [Bbar1 @ H11, Bbar2 + Bbar1 @ H12]
         ])
 
         Cres = np.block([
-            [C1 + D12.dot(T21), D12.dot(T22)],
-            [Dbar21.dot(T11), Cbar2 + Dbar21.dot(T12)]
+            [C1 + D12 @ T21, D12 @ T22],
+            [Dbar21 @ T11, Cbar2 + Dbar21 @ T12]
         ])
 
         Dres = np.block([
-            [D11 + D12.dot(H21), D12.dot(H22)],
-            [Dbar21.dot(H11), Dbar22 + Dbar21.dot(H12)]
+            [D11 + D12 @ H21, D12 @ H22],
+            [Dbar21 @ H11, Dbar22 + Dbar21 @ H12]
         ])
         return StateSpace(Ares, Bres, Cres, Dres, dt)
 
@@ -1381,13 +1381,13 @@ def dynamics(self, t, x, u=None):
         if np.size(x) != self.nstates:
             raise ValueError("len(x) must be equal to number of states")
         if u is None:
-            return self.A.dot(x).reshape((-1,))  # return as row vector
+            return (self.A @ x).reshape((-1,))  # return as row vector
         else:  # received t, x, and u, ignore t
             u = np.reshape(u, (-1, 1))  # force to column in case matrix
             if np.size(u) != self.ninputs:
                 raise ValueError("len(u) must be equal to number of inputs")
-            return self.A.dot(x).reshape((-1,)) \
-                + self.B.dot(u).reshape((-1,))  # return as row vector
+            return (self.A @ x).reshape((-1,)) \
+                + (self.B @ u).reshape((-1,))  # return as row vector
 
     def output(self, t, x, u=None):
         """Compute the output of the system
@@ -1424,13 +1424,13 @@ def output(self, t, x, u=None):
             raise ValueError("len(x) must be equal to number of states")
 
         if u is None:
-            return self.C.dot(x).reshape((-1,))  # return as row vector
+            return (self.C @ x).reshape((-1,))  # return as row vector
         else:  # received t, x, and u, ignore t
             u = np.reshape(u, (-1, 1))  # force to a column in case matrix
             if np.size(u) != self.ninputs:
                 raise ValueError("len(u) must be equal to number of inputs")
-            return self.C.dot(x).reshape((-1,)) \
-                + self.D.dot(u).reshape((-1,))  # return as row vector
+            return (self.C @ x).reshape((-1,)) \
+                + (self.D @ u).reshape((-1,))  # return as row vector
 
     def _isstatic(self):
         """True if and only if the system has no dynamics, that is,
@@ -1623,7 +1623,7 @@ def _rss_generate(states, inputs, outputs, cdtype, strictly_proper=False):
     while True:
         T = randn(states, states)
         try:
-            A = dot(solve(T, A), T)  # A = T \ A * T
+            A = solve(T, A) @ T  # A = T \ A @ T
             break
         except LinAlgError:
             # In the unlikely event that T is rank-deficient, iterate again.

control/timeresp.py

@@ -1972,7 +1972,7 @@ def _ideal_tfinal_and_dt(sys, is_step=True):
         # Incorporate balancing to outer factors
         l[perm, :] *= np.reciprocal(sca)[:, None]
         r[perm, :] *= sca[:, None]
-        w, v = sys_ss.C.dot(r), l.T.conj().dot(sys_ss.B)
+        w, v = sys_ss.C @ r, l.T.conj() @ sys_ss.B
 
         origin = False
         # Computing the "size" of the response of each simple mode