Merge pull request #111 from dynamicslab/control-fix

[Bug fix] Multivariable control inputs for discrete systems

pysindy/optimizers/constrained_sr3.py

@@ -20,7 +20,7 @@ class ConstrainedSR3(SR3):
         0.5\\|y-Xw\\|^2_2 + \\lambda \\times R(v)
         + (0.5 / nu)\\|w-v\\|^2_2
 
-        subject to
+        \\text{subject to}
 
     .. math::
 
@@ -84,7 +84,6 @@ class ConstrainedSR3(SR3):
         "feature" indicates that the constraints are grouped by library feature:
         the first ``n_targets`` columns correspond to the first library feature,
         the next ``n_targets`` columns to the second library feature, and so on.
-        ""
 
     normalize : boolean, optional (default False)
         This parameter is ignored when fit_intercept is set to False. If True,

pysindy/pysindy.py

@@ -755,7 +755,7 @@ def check_stop_condition(xi):
                         return x[: i + 1]
             else:
                 for i in range(1, t):
-                    x[i] = self.predict(x[i - 1 : i], u=u[i - 1])
+                    x[i] = self.predict(x[i - 1 : i], u=u[i - 1, newaxis])
                     if check_stop_condition(x[i]):
                         return x[: i + 1]
             return x

test/conftest.py

@@ -224,6 +224,25 @@ def logistic_map(x, mu, ui):
     return x, u
 
 
+@pytest.fixture
+def data_discrete_time_c_multivariable():
+    def logistic_map(x, mu, u):
+        return mu * x * (1 - x) + u[0] * u[1]
+
+    n_steps = 100
+    mu = 3.6
+
+    u1 = 0.1 * np.random.randn(n_steps)
+    u2 = 0.1 * np.random.randn(n_steps)
+    u = np.column_stack((u1, u2))
+    x = np.zeros((n_steps))
+    x[0] = 0.5
+    for i in range(1, n_steps):
+        x[i] = logistic_map(x[i - 1], mu, u[i - 1])
+
+    return x, u
+
+
 @pytest.fixture
 def data_discrete_time_multiple_trajectories_c():
     def logistic_map(x, mu, ui):

test/test_sindyc.py

@@ -310,8 +310,15 @@ def test_score_multiple_trajectories(data_multiple_trajctories):
     assert s <= 1
 
 
-def test_fit_discrete_time(data_discrete_time_c):
-    x, u = data_discrete_time_c
+@pytest.mark.parametrize(
+    "data",
+    [
+        pytest.lazy_fixture("data_discrete_time_c"),
+        pytest.lazy_fixture("data_discrete_time_c_multivariable"),
+    ],
+)
+def test_fit_discrete_time(data):
+    x, u = data
 
     model = SINDy(discrete_time=True)
     model.fit(x, u=u)
@@ -322,8 +329,15 @@ def test_fit_discrete_time(data_discrete_time_c):
     check_is_fitted(model)
 
 
-def test_simulate_discrete_time(data_discrete_time_c):
-    x, u = data_discrete_time_c
+@pytest.mark.parametrize(
+    "data",
+    [
+        pytest.lazy_fixture("data_discrete_time_c"),
+        pytest.lazy_fixture("data_discrete_time_c_multivariable"),
+    ],
+)
+def test_simulate_discrete_time(data):
+    x, u = data
     model = SINDy(discrete_time=True)
     model.fit(x, u=u)
     n_steps = x.shape[0]
@@ -334,15 +348,29 @@ def test_simulate_discrete_time(data_discrete_time_c):
     # TODO: implement test using the stop_condition option
 
 
-def test_predict_discrete_time(data_discrete_time_c):
-    x, u = data_discrete_time_c
+@pytest.mark.parametrize(
+    "data",
+    [
+        pytest.lazy_fixture("data_discrete_time_c"),
+        pytest.lazy_fixture("data_discrete_time_c_multivariable"),
+    ],
+)
+def test_predict_discrete_time(data):
+    x, u = data
     model = SINDy(discrete_time=True)
     model.fit(x, u=u)
     assert len(model.predict(x, u=u)) == len(x)
 
 
-def test_score_discrete_time(data_discrete_time_c):
-    x, u = data_discrete_time_c
+@pytest.mark.parametrize(
+    "data",
+    [
+        pytest.lazy_fixture("data_discrete_time_c"),
+        pytest.lazy_fixture("data_discrete_time_c_multivariable"),
+    ],
+)
+def test_score_discrete_time(data):
+    x, u = data
     model = SINDy(discrete_time=True)
     model.fit(x, u=u)
     assert model.score(x, u=u) > 0.75