Merge pull request #1870 from nmegiddo/patch-5

Update utils.py

Author: beat-buesser

art/utils.py

@@ -376,13 +376,13 @@ def projection_l1_1(values: np.ndarray, eps: Union[int, float, np.ndarray]) -> n
     mat = np.zeros((m, 2))
 
     #   if  a_sorted[i, n-1]  >= a_sorted[i, n-2] + eps,  then the projection is  [0,...,0,eps]
-    done = False
+    done = early_done = False
     active = np.array([1] * m)
     after_vec = np.zeros((m, n))
     proj = a_sorted.copy()
     j = n - 2
     while j >= 0:
-        mat[:, 0] = mat[:, 0] + a_sorted[:, j + 1]  # =  sum(a_sorted[: i] :  i = j + 1,...,n-1
+        mat[:, 0] += a_sorted[:, j + 1]  # =  sum(a_sorted[: i] :  i = j + 1,...,n-1
         mat[:, 1] = a_sorted[:, j] * (n - j - 1) + eps
         #  Find the max in each problem  max{ sum{a_sorted[:, i] : i=j+1,..,n-1} , a_sorted[:, j] * (n-j-1) + eps }
         row_maxes = np.max(mat, axis=1)
@@ -396,21 +396,29 @@ def projection_l1_1(values: np.ndarray, eps: Union[int, float, np.ndarray]) -> n
         #  has to be reduced is  delta
         delta = (mat[:, 0] - eps) / (n - j - 1)
         #    The vector of reductions
-        delta_vec = np.array([delta] * (n - j - 1))
-        delta_vec = np.transpose(delta_vec)
+        delta_vec = np.transpose(np.array([delta] * (n - j - 1)))
         #   The sub-vectors:  a_sorted[:, (j+1):]
         a_sub = a_sorted[:, (j + 1) :]
         #   After reduction by delta_vec
         a_after = a_sub - delta_vec
         after_vec[:, (j + 1) :] = a_after
-        proj = (act_multiplier * after_vec) + ((1 - act_multiplier) * proj)
+        proj += act_multiplier * (after_vec - proj)
         active = active * ind_set
         if sum(active) == 0:
-            done = True
+            done = early_done = True
             break
         j -= 1
+    if not early_done:
+        delta = (mat[:, 0] + a_sorted[:, 0] - eps) / n
+        ind_set = np.sign(np.maximum(delta, 0))
+        act_multiplier = ind_set * active
+        act_multiplier = np.transpose([np.transpose(act_multiplier)] * n)
+        delta_vec = np.transpose(np.array([delta] * n))
+        a_after = a_sorted - delta_vec
+        proj += act_multiplier * (a_after - proj)
+        done = True
     if not done:
-        proj = active * a_sorted + (1 - active) * proj
+        proj = active * (a_sorted - proj)
 
     for i in range(m):
         proj[i, :] = proj[i, a_argsort_inv[i, :]]
@@ -461,7 +469,7 @@ def projection_l1_2(values: np.ndarray, eps: Union[int, float, np.ndarray]) -> n
         mat0[:, 1] = np.min(mat, axis=1)
         min_t = np.max(mat0, axis=1)
         if np.max(min_t) < 1e-8:
-            break
+            continue
         row_sums = row_sums - a_var[:, j] * (n - j)
         a_var[:, (j + 1) :] = a_var[:, (j + 1) :] - np.matmul(min_t.reshape((m, 1)), np.ones((1, n - j - 1)))
         a_var[:, j] = a_var[:, j] - min_t