Simplex Algorithm
11-12-2023, 05:11 PM (This post was last modified: 11-12-2023 11:23 PM by Albert Chan.)
Post: #6
 Albert Chan Senior Member Posts: 2,636 Joined: Jul 2018
RE: Simplex Algorithm
(11-12-2023 12:38 AM)Albert Chan Wrote:  Fastest is remove excess variables, turning it into slack variables.

Note that we scale pivot as 1 before clearing other cells.
This does not affect other constraints, not even cost functions.
In other words, row operations are safe.

Time (including pivot operations) = 0.0075 s
Production code need pivots selection, but should not cost much.

I don't know a good name for this algorithm, so I just use remove_eq_c(a, n)
It remove top n "=" constraints. (n variables, if constraints are independent)
In other words, only "≤" constraints remain.

simplex_le(a, ±n) ≡ simplex_le(remove_eq_c(a,n), ±∞)

simplex_le(a, dir) only had 1 change. Instead of more constraints, it now do less variables.

< IF b<len(a) THEN a := extend(-a[1 .. b],a) END
> IF b<len(a) THEN a := remove_eq_c(a, b) END

Code:
pivot_index(a) := BEGIN         LOCAL b := contains((a:=abs(a)), 1.0);     IF b THEN RETURN b END     IF (b:=max(a)) THEN b := contains(a,b) END     RETURN b; END; remove_eq_c(a, n) := BEGIN   LOCAL p, r, c := len(a[1]);   FOR r FROM n DOWNTO 1 DO     p := pivot_index(suppress(a[r],c));     IF p==0 and a[r][c] THEN error(a) END     a := p ? pivot((a[r]/=a[r][p]),r,p) : delrows(a,r);   END   return a; END; simplex_le(a, dir) := BEGIN  LOCAL b := abs(dir);  IF b<len(a) THEN a := remove_eq_c(a, b) END  b := dim(a);  a := transpose(extend(     [col(a,1 .. b(2)-1)],     [col(identity(b(1)),1 .. b(1)-1)],     [col(a,b(2))]  ));  IF dir>0 THEN a[0] := -a[0] END /* maximize */  a := simplex2(a,range(b(2),sum(b)-1),sum(b)-2,0,0);  IF dir>0 THEN a[1] := -a[1] END /* undo flip */  return a; END;

Benchmark example:

> m := [[3,2,1,10], [2,5,3,15], [-2,-3,-4,0]]
> time(r := simplex_le(m, -2))                     → 0.0077
> r

$$[\displaystyle \frac{-130}{7},[\frac{15}{7},0,\frac{25}{7},\frac{15}{7},0],\left(\begin{array}{cccccc} 1 & \frac{1}{7} & 0 & 1 & \frac{1}{7} & \frac{15}{7} \\ 0 & \frac{11}{7} & 1 & 0 & \frac{11}{7} & \frac{25}{7} \\ 0 & \frac{25}{7} & 0 & 0 & \frac{25}{7} & \frac{130}{7} \end{array}\right) ]$$

Update Nov 12, 2023: simplex_le directly call simplex2, instead of wrapper simplex.