Fix error when removing unused variables

Author: jsiirola

pyomo/repn/plugins/standard_form.py

@@ -472,24 +472,22 @@ def write(self, model):
         # at the index pointer list (an O(num_var) operation).
         c_ip = c.indptr
         A_ip = A.indptr
-        active_var_idx = list(
-            filter(
-                lambda i: A_ip[i] != A_ip[i + 1] or c_ip[i] != c_ip[i + 1],
-                range(len(columns)),
-            )
-        )
-        nCol = len(active_var_idx)
+        active_var_mask = (A_ip[1:] > A_ip[:-1]) | (c_ip[1:] > c_ip[:-1])
+
+        # Masks on NumPy arrays are very fast.  Build the reduced A
+        # indptr and then check if we actually have to manipulate the
+        # columns
+        augmented_mask = np.concatenate((active_var_mask, [True]))
+        reduced_A_indptr = A.indptr[augmented_mask]
+        nCol = len(reduced_A_indptr) - 1
         if nCol != len(columns):
-            # Note that the indptr can't just use range() because a var
-            # may only appear in the objectives or the constraints.
-            columns = list(map(columns.__getitem__, active_var_idx))
-            active_var_idx.append(c.indptr[-1])
+            columns = [v for k, v in zip(active_var_mask, columns) if k]
             c = scipy.sparse.csc_array(
-                (c.data, c.indices, c.indptr.take(active_var_idx)), [c.shape[0], nCol]
+                (c.data, c.indices, c.indptr[augmented_mask]), [c.shape[0], nCol]
             )
-            active_var_idx[-1] = A.indptr[-1]
+            # active_var_idx[-1] = len(columns)
             A = scipy.sparse.csc_array(
-                (A.data, A.indices, A.indptr.take(active_var_idx)), [A.shape[0], nCol]
+                (A.data, A.indices, reduced_A_indptr), [A.shape[0], nCol]
             )
 
         if self.config.nonnegative_vars:

pyomo/repn/tests/test_standard_form.py

@@ -43,6 +43,19 @@ def test_linear_model(self):
         self.assertTrue(np.all(repn.A == np.array([[-1, -2, 0], [0, 1, 4]])))
         self.assertTrue(np.all(repn.rhs == np.array([-3, 5])))
 
+    def test_almost_dense_linear_model(self):
+        m = pyo.ConcreteModel()
+        m.x = pyo.Var()
+        m.y = pyo.Var([1, 2, 3])
+        m.c = pyo.Constraint(expr=m.x + 2 * m.y[1] + 4 * m.y[3] >= 10)
+        m.d = pyo.Constraint(expr=5 * m.x + 6 * m.y[1] + 8 * m.y[3] <= 20)
+
+        repn = LinearStandardFormCompiler().write(m)
+
+        self.assertTrue(np.all(repn.c == np.array([0, 0, 0])))
+        self.assertTrue(np.all(repn.A == np.array([[-1, -2, -4], [5, 6, 8]])))
+        self.assertTrue(np.all(repn.rhs == np.array([-10, 20])))
+
     def test_linear_model_row_col_order(self):
         m = pyo.ConcreteModel()
         m.x = pyo.Var()