improve tests by scripting julia code to mimic python tests, and seri…

…alize with json file. Fix scalar multiplication
MythasNauveili · Mar 26, 2025 · Mar 26, 2025 · Mar 26, 2025 · Mar 28, 2025 · Mar 29, 2025
commit 30c05b93e957442ac85f3a26a9207bf12bfe7690
diff --git a/control/dde.py b/control/dde.py
@@ -244,6 +244,12 @@ def f(t, x):
     return xout, yout
 
 
+
+def dde23_solver(delay_sys, dt, T, U, X0, xout, yout):
+    # Plans for implementing a python version matlab dde23 solver
+    raise NotImplementedError
+
+
 def linear_interp_u(t, T, U):
     """
     Linearly interpolate the input U at time t.

diff --git a/control/delaylti.py b/control/delaylti.py
@@ -70,6 +70,22 @@ class DelayLTI(LTI):
         Create a DelayLTI system from a TransferFunction system.
     size()
         Return the number of outputs and inputs.
+
+    Notes
+    -----
+    The preferred way to create a DelayLTI object is by multiplying a transfert
+    function object with a delay(tau) or exp(-tau*s). For example
+
+    >>> ct.tf([1], [1,1]) * delay(1.5)
+
+    Or
+
+    >>> s = ct.tf('s')
+    >>> ct.tf([1], [1,1]) * exp(-1.5*s)
+
+    It's possible to MIMO delayed systems from arrays of SISO delayed systems,
+    see function mimo_delay
+
     """
 
     def __init__(self, P: PartitionedStateSpace, tau = None, **kwargs):
@@ -154,6 +170,10 @@ def size(self):
         """Return the number of outputs and inputs."""
         return (self.noutputs, self.ninputs)
 
+    # Poles and zeros functions for DelayLTI are not supported by julia ControlSystems.jl
+    # might not be accurate for delayLTI, to be discussed
+    # Here the poles and zeros computed are the ones of the system without taking the delay
+
     def poles(self):
         """Compute the poles of a delay lti system."""
 
@@ -229,9 +249,10 @@ def __mul__(self, other):
         """
 
         if isinstance(other, (int, float, complex)):
-            new_C = np.block([[self.P.C1 * other], [self.P.C2]])
-            new_D = np.block([[self.P.D11 * other, self.P.D12 * other], [self.P.D21, self.P.D22]])
-            new_P = PartitionedStateSpace(ss(self.P.A, self.P.B, new_C, new_D), self.P.nu1, self.P.ny1)
+            new_B = np.hstack([self.P.B1 * other, self.P.B2])
+            new_D = np.block([[self.P.D11 * other, self.P.D12], [self.P.D21 * other, self.P.D22]])
+
+            new_P = PartitionedStateSpace(ss(self.P.A, new_B, self.P.C, new_D), self.P.nu1, self.P.ny1)
             return DelayLTI(new_P, self.tau)
 
         elif isinstance(other, DelayLTI):
@@ -249,10 +270,14 @@ def __mul__(self, other):
         else:
             raise TypeError("Unsupported operand type(s) for *: '{}' and '{}'".format(type(self), type(other)))
 
+
     def __rmul__(self, other):
-
         if isinstance(other, (int, float, complex)):
-            return self * other
+            new_C = np.hstack([self.P.C1 * other, self.P.C2])
+            new_D = np.block([[self.P.D11 * other, self.P.D12 * other], [self.P.D21, self.P.D22]])
+
+            new_P = PartitionedStateSpace(ss(self.P.A, self.P.B, new_C, new_D), self.P.nu1, self.P.ny1)
+            return DelayLTI(new_P, self.tau)
 
         elif isinstance(other, TransferFunction):
             dlti = tf2dlti(other)
@@ -467,6 +492,14 @@ def __str__(self):
     def __repr__(self):
          return (f"{type(self).__name__}("
                  f"P={self.P.__repr__()}, tau={self.tau.__repr__()})")
+
+    # Some functions from LTI that are not currently supported with delayLTI
+    def to_ss(self, *args, **kwargs):
+        return NotImplementedError
+
+    def to_tf(self, *args, **kwargs):
+        raise NotImplementedError
+
 
 
 def delay(tau):

diff --git a/control/julia_utils.py b/control/julia_utils.py
@@ -0,0 +1,55 @@
+import numpy as np
+import json
+
+from pathlib import Path
+
+def _recursive_reshape(node):
+    if isinstance(node, dict) and ("data" in node) and ("dim" in node):
+        data = node["data"]
+        dim = tuple(node["dim"])
+
+        array_data = np.array(data)
+
+        if len(dim) == 1:
+            return array_data
+
+        elif len(dim) > 1 and (np.shape(array_data) == dim) and (dim[0] != dim[1]):
+            return array_data
+        else:
+            return np.transpose(array_data)
+
+    elif isinstance(node, dict):
+        new_node = {}
+        for key, value in node.items():
+            new_node[key] = _recursive_reshape(value)
+        return new_node
+
+    elif isinstance(node, list) and any(isinstance(item, dict) for item in node):
+        new_node = []
+        for i, item in enumerate(node):
+            new_node[i] = _recursive_reshape(item)
+        return new_node
+
+    else:
+        return node
+
+
+def load_julia_results(json_path):
+    with open(json_path, "r") as f:
+        json_content = json.load(f)
+
+    reshaped_julia = _recursive_reshape(json_content)
+    return reshaped_julia
+
+
+def assert_delayLTI(dlti, julia_results):
+    assert np.allclose(dlti.P.A, julia_results["A"])
+    assert np.allclose(dlti.P.B, julia_results["B"])
+    assert np.allclose(dlti.P.C, julia_results["C"])
+    assert np.allclose(dlti.P.D, julia_results["D"])
+    assert np.allclose(dlti.tau, julia_results["tau"])
+
+
+# Load julia results file
+script_dir = Path(__file__).parent
+julia_json = load_julia_results(f"{script_dir}/tests/julia_results.json")