neil-lindquist/linear-programming

A Common Lisp library for solving linear programming problems

Common Lisp Linear Programming

This is a Common Lisp library for solving linear programming problems. It's designed to provide a high-level and ergonomic API for specifying linear programming problems as lisp expressions.

The core library is implemented purely in Common Lisp with only a few community-standard libraries as dependencies (ASDF, Alexandria, Iterate). However, the solver is designed to support alternative backends without any change to the user's code. Currently, there is a backend for the GNU Linear Programming Kit (GLPK).

Installation

The linear-programming library is avalible in both the main Quicklisp distribution and Ultralisp, so it can loaded with with (ql:quickload :linear-programming). You can check that it works by running (asdf:test-system :linear-programming).

If you are not using Quicklisp, place this repository, Alexandria, and Iterate somewhere where ASDF can find them. Then, it can be loaded with (asdf:load-system :linear-programming) and tested as above.

Usage

See neil-lindquist.github.io/linear-programming/ for further documentation.

Consider the following linear programming problem.

maximize x + 4y + 3z
such that

• 2x + y ≤ 8
• y + z ≤ 7
• x, y, z ≥ 0

First, the problem needs to be specified. Problems are specified with a simple DSL, as described in the syntax reference.

(use-package :linear-programming)

(defvar problem (parse-linear-problem '(max (= w (+ x (* 4 y) (* 3 z))))
                                      '((<= (+ (* 2 x) y) 8)
                                        (<= (+ y z) 7))))

Once the problem is created, it can be solved with the simplex method.

(defvar solution (solve-problem problem))

Finally, the optimal tableau can be inspected to get the resulting objective function, decision variables, and reduced-costs (i.e. the shadow prices for the variable's lower bounds).

(format t "Objective value solution: ~A~%" (solution-variable solution 'w))
(format t "x = ~A (reduced cost: ~A)~%" (solution-variable solution 'x) (solution-reduced-cost solution 'x))
(format t "y = ~A (reduced cost: ~A)~%" (solution-variable solution 'y) (solution-reduced-cost solution 'y))
(format t "z = ~A (reduced cost: ~A)~%" (solution-variable solution 'z) (solution-reduced-cost solution 'z))

;; ==>
;; Objective value solution: 57/2
;; x = 1/2 (reduced cost: 0)
;; y = 7 (reduced cost: 0)
;; z = 0 (reduced cost: 1/2)

Alternatively, the with-solution-variables and with-solved-problem macros simplify some steps and binds the solution variables in their bodies.

(with-solution-variables (w x y z) solution
  (format t "Objective value solution: ~A~%" w)
  (format t "x = ~A (reduced cost: ~A)~%" x (reduced-cost x))
  (format t "y = ~A (reduced cost: ~A)~%" y (reduced-cost y))
  (format t "z = ~A (reduced cost: ~A)~%" z (reduced-cost z)))

;; ==>
;; Objective value solution: 57/2
;; x = 1/2 (reduced cost: 0)
;; y = 7 (reduced cost: 0)
;; z = 0 (reduced cost: 1/2)


(with-solved-problem ((max (= w (+ x (* 4 y) (* 3 z))))
                      (<= (+ (* 2 x) y) 8)
                      (<= (+ y z) 7))
  (format t "Objective value solution: ~A~%" w)
  (format t "x = ~A (reduced cost: ~A)~%" x (reduced-cost x))
  (format t "y = ~A (reduced cost: ~A)~%" y (reduced-cost y))
  (format t "z = ~A (reduced cost: ~A)~%" z (reduced-cost z)))

;; ==>
;; Objective value solution: 57/2
;; x = 1/2 (reduced cost: 0)
;; y = 7 (reduced cost: 0)
;; z = 0 (reduced cost: 1/2)

linear-programming's Issues

Feasible integer problem reported as unbounded

Thank you for fixing #7 , but now a slightly larger problem fails (for a different reason):

(lp:with-solved-problem 
	     ((MAX (= W (+ (* 1 FUEL))))
	      (INTEGER PSHF T65 FUEL T63 HKGWZ T66 GPVTF T68 QDVJ T64 KHKGT T69 XJWVT T67
		       DCFZ T62 NZVS T61 ORE)
	      (LINEAR-PROGRAMMING/PROBLEM:BOUNDS (0 ORE 1000000000000)) (= (+ (* -1 PSHF) (* 7 T65)) 0)
	      (= (+ (* -1 FUEL) (* 1 T63)) 0) (= (+ (* -1 HKGWZ) (* 5 T66)) 0)
	      (= (+ (* -1 GPVTF) (* 2 T68)) 0) (= (+ (* -1 QDVJ) (* 9 T64)) 0)
	      (= (+ (* -1 KHKGT) (* 8 T69)) 0) (= (+ (* -1 XJWVT) (* 2 T67)) 0)
	      (= (+ (* -1 DCFZ) (* 6 T62)) 0) (= (+ (* -1 NZVS) (* 5 T61)) 0)
	      (<=
	       (+ (* -1 PSHF) (* 8 T17) (* 7 T20) (* 10 T22) (* 8 T26) (* 7 T29) (* 10 T31)
		  (* 8 T64) (* 7 T67) (* 10 T69))
	       0)
	      (<=
	       (+ (* -1 HKGWZ) (* 48 T16) (* 12 T17) (* 5 T22) (* 48 T25) (* 12 T26)
		  (* 5 T31) (* 48 T63) (* 12 T64) (* 5 T69))
	       0)
	      (<=
	       (+ (* -1 GPVTF) (* 9 T16) (* 1 T17) (* 9 T25) (* 1 T26) (* 9 T63) (* 1 T64))
	       0)
	      (<= (+ (* -1 QDVJ) (* 1 T16) (* 1 T25) (* 1 T63)) 0)
	      (<= (+ (* -1 KHKGT) (* 5 T16) (* 5 T25) (* 5 T63)) 0)
	      (<= (+ (* -1 XJWVT) (* 44 T16) (* 44 T25) (* 44 T63)) 0)
	      (<=
	       (+ (* -1 DCFZ) (* 7 T20) (* 3 T22) (* 7 T29) (* 3 T31) (* 7 T67) (* 3 T69))
	       0)
	      (<=
	       (+ (* -1 NZVS) (* 29 T16) (* 7 T22) (* 29 T25) (* 7 T31) (* 29 T63)
		  (* 7 T69))
	       0)
	      (<=
	       (+ (* -1 ORE) (* 10 T1) (* 1 T2) (* 9 T7) (* 8 T8) (* 7 T9) (* 157 T14)
		  (* 165 T15) (* 179 T18) (* 177 T19) (* 165 T21) (* 157 T23) (* 165 T24)
		  (* 179 T27) (* 177 T28) (* 165 T30) (* 139 T36) (* 144 T37) (* 145 T40)
		  (* 176 T43) (* 171 T44) (* 114 T46) (* 189 T53) (* 121 T58) (* 157 T61)
		  (* 165 T62) (* 179 T65) (* 177 T66) (* 165 T68))
	       0))
	   w)

This raises an UNBOUNDED-PROBLEM-ERROR but only with the default Simplex backend. GLPK backend works as intended (w = 82892753).

Feasible problem reported unbounded

Using the code from Ultralisp, the following raises an unbounded problem condition:

(with-input-from-string (in "((MIN (+ (* 0.6861807 A) (* 1 B))) (>= (+ (* 1 B) (* 0.6861807 A)) 0.9372585) (>= (+ (* 1 B) (* 0.7776901 A)) 0.7461006) (>= (+ (* 1 B) (* 0.14247864 A)) 0.38555977))") (lp:solve-problem (lp:read-sexp in)))

That should not be possible (due to the way the problem is constructed). Here's an optimal solution found by LPSolve that I have checked is correct.

a=1.01470781626246 b=0.240985580341555

Am I using the API wrong or is this a bug?

Supporting GLPK

Hi, the README mentions

If there is interest in a high performance backend, please open an issue; the library is setup to support a high performance backend, but there has yet to be any interest in it.

and I'd like to say that I'd be interested in it, and helping to get it done. A while ago, I was playing with cl-glpk, and got it working to some degree, so I suppose that implementing a GLPK backend for this library would go along the same lines.

Additional Examples

Nice system, thanks. Perhaps some additional documentation and examples, especially examples, would be helpful. For example one demonstrating how to do linear regression?

Array index-out-of-bounds when solving ILP

Using the latest code from Ultralisp:

(lp:with-solved-problem 
	     ((MIN (= W (+ (* 1 ORE))))
	      (INTEGER FUEL T185 E T184 D T183 C T182 B T181 A T180 ORE)
	      (LINEAR-PROGRAMMING/PROBLEM:BOUNDS (1 FUEL 1))
	      (= (+ (* -1 FUEL) (* 1 T185)) 0) (= (+ (* -1 E) (* 1 T184)) 0)
	      (= (+ (* -1 D) (* 1 T183)) 0) (= (+ (* -1 C) (* 1 T182)) 0)
	      (= (+ (* -1 B) (* 1 T181)) 0) (= (+ (* -1 A) (* 10 T180)) 0)
	      (<= (+ (* -1 E) (* 1 T185)) 0) (<= (+ (* -1 D) (* 1 T184)) 0)
	      (<= (+ (* -1 C) (* 1 T183)) 0) (<= (+ (* -1 B) (* 1 T182)) 0)
	      (<= (+ (* -1 A) (* 7 T182) (* 7 T183) (* 7 T184) (* 7 T185)) 0)
	      (<=
	       (+ (* -1 ORE) (* 171 T1) (* 114 T3) (* 189 T10) (* 121 T15) (* 156 T18)
		  (* 185 T52) (* 111 T54) (* 141 T63) (* 156 T72) (* 185 T106) (* 111 T108)
		  (* 141 T117) (* 156 T126) (* 185 T160) (* 111 T162) (* 141 T171)
		  (* 10 T180) (* 1 T181))
	       0))
	   w)

fails with the error "Invalid index 37 for (SIMPLE-ARRAY T (14 37)), should be a non-negative integer below 37."
If I use the GLPK backend, I get the expected solution of w = 31. So something must be wrong in your Simplex implementation, looks like it might just be an off-by-one :)

Solution does not appear to be optimal

Hi,

Thanks so much for fixing issue #3!

After running the solver on the Stigler Diet problem, I noticed a few issues.

(defparameter data
  '(("Wheat Flour (Enriched)" "10 lb." 36 44.7 1411 2 365 0 55.4 33.3 441 0)
    ("Macaroni" "1 lb." 14.1 11.6 418 0.7 54 0 3.2 1.9 68 0)
    ("Wheat Cereal (Enriched)" "28 oz." 24.2 11.8 377 14.4 175 0 14.4 8.8 114 0)
    ("Corn Flakes" "8 oz." 7.1 11.4 252 0.1 56 0 13.5 2.3 68 0)
    ("Corn Meal" "1 lb." 4.6 36.0 897 1.7 99 30.9 17.4 7.9 106 0)
    ("Hominy Grits" "24 oz." 8.5 28.6 680 0.8 80 0 10.6 1.6 110 0)
    ("Rice" "1 lb." 7.5 21.2 460 0.6 41 0 2 4.8 60 0)
    ("Rolled Oats" "1 lb." 7.1 25.3 907 5.1 341 0 37.1 8.9 64 0)
    ("White Bread (Enriched)" "1 lb." 7.9 15.0 488 2.5 115 0 13.8 8.5 126 0)
    ("Whole Wheat Bread" "1 lb." 9.1 12.2 484 2.7 125 0 13.9 6.4 160 0)
    ("Rye Bread" "1 lb." 9.1 12.4 439 1.1 82 0 9.9 3 66 0)
    ("Pound Cake" "1 lb." 24.8 8.0 130 0.4 31 18.9 2.8 3 17 0)
    ("Soda Crackers" "1 lb." 15.1 12.5 288 0.5 50 0 0 0 0 0)
    ("Milk" "1 qt." 11 6.1 310 10.5 18 16.8 4 16 7 177)
    ("Evaporated Milk (can)" "14.5 oz." 6.7 8.4 422 15.1 9 26 3 23.5 11 60)
    ("Butter" "1 lb." 30.8 10.8 9 0.2 3 44.2 0 0.2 2 0)
    ("Oleomargarine" "1 lb." 16.1 20.6 17 0.6 6 55.8 0.2 0 0 0)
    ("Eggs" "1 doz." 32.6 2.9 238 1.0 52 18.6 2.8 6.5 1 0)
    ("Cheese (Cheddar)" "1 lb." 24.2 7.4 448 16.4 19 28.1 0.8 10.3 4 0)
    ("Cream" "1/2 pt." 14.1 3.5 49 1.7 3 16.9 0.6 2.5 0 17)
    ("Peanut Butter" "1 lb." 17.9 15.7 661 1.0 48 0 9.6 8.1 471 0)
    ("Mayonnaise" "1/2 pt." 16.7 8.6 18 0.2 8 2.7 0.4 0.5 0 0)
    ("Crisco" "1 lb." 20.3 20.1 0 0 0 0 0 0 0 0)
    ("Lard" "1 lb." 9.8 41.7 0 0 0 0.2 0 0.5 5 0)
    ("Sirloin Steak" "1 lb." 39.6 2.9 166 0.1 34 0.2 2.1 2.9 69 0)
    ("Round Steak" "1 lb." 36.4 2.2 214 0.1 32 0.4 2.5 2.4 87 0)
    ("Rib Roast" "1 lb." 29.2 3.4 213 0.1 33 0 0 2 0 0)
    ("Chuck Roast" "1 lb." 22.6 3.6 309 0.2 46 0.4 1 4 120 0)
    ("Plate" "1 lb." 14.6 8.5 404 0.2 62 0 0.9 0 0 0)
    ("Liver (Beef)" "1 lb." 26.8 2.2 333 0.2 139 169.2 6.4 50.8 316 525)
    ("Leg of Lamb" "1 lb." 27.6 3.1 245 0.1 20 0 2.8 3.9 86 0)
    ("Lamb Chops (Rib)" "1 lb." 36.6 3.3 140 0.1 15 0 1.7 2.7 54 0)
    ("Pork Chops" "1 lb." 30.7 3.5 196 0.2 30 0 17.4 2.7 60 0)
    ("Pork Loin Roast" "1 lb." 24.2 4.4 249 0.3 37 0 18.2 3.6 79 0)
    ("Bacon" "1 lb." 25.6 10.4 152 0.2 23 0 1.8 1.8 71 0)
    ("Ham smoked" "1 lb." 27.4 6.7 212 0.2 31 0 9.9 3.3 50 0)
    ("Salt Pork" "1 lb." 16 18.8 164 0.1 26 0 1.4 1.8 0 0)
    ("Roasting Chicken" "1 lb." 30.3 1.8 184 0.1 30 0.1 0.9 1.8 68 46)
    ("Veal Cutlets" "1 lb." 42.3 1.7 156 0.1 24 0 1.4 2.4 57 0)
    ("Salmon Pink (can)" "16 oz." 13 5.8 705 6.8 45 3.5 1 4.9 209 0)
    ("Apples" "1 lb." 4.4 5.8 27 0.5 36 7.3 3.6 2.7 5 544)
    ("Bananas" "1 lb." 6.1 4.9 60 0.4 30 17.4 2.5 3.5 28 498)
    ("Lemons" "1 doz." 26 1.0 21 0.5 14 0 0.5 0 4 952)
    ("Oranges" "1 doz." 30.9 2.2 40 1.1 18 11.1 3.6 1.3 10 1998)
    ("Green Beans" "1 lb." 7.1 2.4 138 3.7 80 69 4.3 5.8 37 862)
    ("Cabbage" "1 lb." 3.7 2.6 125 4.0 36 7.2 9 4.5 26 5369)
    ("Carrots" "1 bunch" 4.7 2.7 73 2.8 43 188.5 6.1 4.3 89 608)
    ("Celery" "1 stalk" 7.3 0.9 51 3.0 23 0.9 1.4 1.4 9 313)
    ("Lettuce" "1 head" 8.2 0.4 27 1.1 22 112.4 1.8 3.4 11 449)
    ("Onions" "1 lb." 3.6 5.8 166 3.8 59 16.6 4.7 5.9 21 1184)
    ("Potatoes" "15 lb." 34 14.3 336 1.8 118 6.7 29.4 7.1 198 2522)
    ("Spinach" "1 lb." 8.1 1.1 106 0 138 918.4 5.7 13.8 33 2755)
    ("Sweet Potatoes" "1 lb." 5.1 9.6 138 2.7 54 290.7 8.4 5.4 83 1912)
    ("Peaches (can)" "No. 2 1/2" 16.8 3.7 20 0.4 10 21.5 0.5 1 31 196)
    ("Pears (can)" "No. 2 1/2" 20.4 3.0 8 0.3 8 0.8 0.8 0.8 5 81)
    ("Pineapple (can)" "No. 2 1/2" 21.3 2.4 16 0.4 8 2 2.8 0.8 7 399)
    ("Asparagus (can)" "No. 2" 27.7 0.4 33 0.3 12 16.3 1.4 2.1 17 272)
    ("Green Beans (can)" "No. 2" 10 1.0 54 2 65 53.9 1.6 4.3 32 431)
    ("Pork and Beans (can)" "16 oz." 7.1 7.5 364 4 134 3.5 8.3 7.7 56 0)
    ("Corn (can)" "No. 2" 10.4 5.2 136 0.2 16 12 1.6 2.7 42 218)
    ("Peas (can)" "No. 2" 13.8 2.3 136 0.6 45 34.9 4.9 2.5 37 370)
    ("Tomatoes (can)" "No. 2" 8.6 1.3 63 0.7 38 53.2 3.4 2.5 36 1253)
    ("Tomato Soup (can)" "10 1/2 oz." 7.6 1.6 71 0.6 43 57.9 3.5 2.4 67 862)
    ("Peaches Dried" "1 lb." 15.7 8.5 87 1.7 173 86.8 1.2 4.3 55 57)
    ("Prunes Dried" "1 lb." 9 12.8 99 2.5 154 85.7 3.9 4.3 65 257)
    ("Raisins Dried" "15 oz." 9.4 13.5 104 2.5 136 4.5 6.3 1.4 24 136)
    ("Peas Dried" "1 lb." 7.9 20.0 1367 4.2 345 2.9 28.7 18.4 162 0)
    ("Lima Beans Dried" "1 lb." 8.9 17.4 1055 3.7 459 5.1 26.9 38.2 93 0)
    ("Navy Beans Dried" "1 lb." 5.9 26.9 1691 11.4 792 0 38.4 24.6 217 0)
    ("Coffee" "1 lb." 22.4 0 0 0 0 0 4 5.1 50 0)
    ("Tea" "1/4 lb." 17.4 0 0 0 0 0 0 2.3 42 0)
    ("Cocoa" "8 oz." 8.6 8.7 237 3 72 0 2 11.9 40 0)
    ("Chocolate" "8 oz." 16.2 8.0 77 1.3 39 0 0.9 3.4 14 0)
    ("Sugar" "10 lb." 51.7 34.9 0 0 0 0 0 0 0 0)
    ("Corn Syrup" "24 oz." 13.7 14.7 0 0.5 74 0 0 0 5 0)
    ("Molasses" "18 oz." 13.6 9.0 0 10.3 244 0 1.9 7.5 146 0)
    ("Strawberry Preserves" "1 lb." 20.5 6.4 11 0.4 7 0.2 0.2 0.4 3 0)))

(defparameter nutrients
  '(("Calories (1000s)" 3)
    ("Protein (grams)" 70)
    ("Calcium (grams)" 0.8)
    ("Iron (mg)" 12)
    ("Vitamin A (1000 IU)" 5)
    ("Vitamin B1 (mg)" 1.8)
    ("Vitamin B2 (mg)" 2.7)
    ("Niacin (mg)" 18)
    ("Vitamin C (mg)" 75)))

(defparameter objective-function
  `(min (= w (+ ,@(mapcar 'make-symbol (mapcar 'first data))))))

(defparameter constraints
  (let ((constraints (list)))
    (do ((i 0 (1+ i)))
        ((= i (length nutrients)) i)
      (push `(<= ,(second (nth i nutrients))
                 (+ ,@(mapcar (lambda (d) `(* ,(nth (+ i 3) d)
                                              ,(make-symbol (first d)))) data)))
            constraints))
    constraints))

(defparameter problem
  (parse-linear-problem
   objective-function
   constraints))

(defparameter solution
  (solve-problem problem))

(defparameter solution-vars
  (map 'vector
       (lambda (var) (list var (solution-variable solution var)))
       (problem-vars problem)))

(defparameter non-zero-solution-vars
  (remove-if (lambda (x) (zerop (second x))) solution-vars))

(defparameter optimal-annual-price
  (* (loop for x across (map 'vector 'second non-zero-solution-vars) sum x) 365))

(defparameter data-hash
  (let ((hash-table (make-hash-table)))
    (dolist (datum data)
      (setf (gethash (first datum) hash-table) datum))
    hash-table))

(defparameter nutrient-sums
  (let* ((num-nutrients (length nutrients))
         (sums (make-array num-nutrients)))
    (do ((i 0 (1+ i)))
        ((= i num-nutrients) i)
      (setf (elt sums i)
            (loop for var across non-zero-solution-vars sum (* (second var)
                                                               (elt (gethash
                                                                     (symbol-name (first var))
                                                                     data-hash)
                                                                    (+ i 3))))))
    sums))

1. If I inspect problem by entering it in the repl, I notice that :vars contains duplicate variables.

2. If I inspect solution-vars in the repl, I notice that it, too, has duplicate variables.

3. If I inspect non-zero-solution-vars in the repl, I notice that it, too, has duplicate variables (namely #:|Wheat Flour (Enriched)|).

#((#:|Wheat Flour (Enriched)| 0.06711409)
  (#:|Cheese (Cheddar)| 0.04878049)
  (#:|Tomato Soup (can)| 0.08635578)
  (#:|Liver (Beef)| 0.053149607)
  (#:|Cabbage| 0.013969081)
  (#:|Peanut Butter| 0.03821656)
  (#:|Wheat Flour (Enriched)| 0.032490972)
  (#:|Corn Flakes| 0.21428573)
  (#:|Navy Beans Dried| 0.041395623))

4. If I inspect optimal-annual-price, I notice that it is 217.45163, whereas in the output of GLOP (https://developers.google.com/optimization/lp/glop) it is 39.66:

Wheat Flour (Enriched) = 0.029519
Liver (Beef) = 0.001893
Cabbage = 0.011214
Spinach = 0.005008
Navy Beans, Dried = 0.061029
Optimal annual price: $39.66

5. If I inspect nutrient-sums, I see the following:

#(9.26114 337.23376 1.6490852 95.506386 15.464223 11.174599 8.610146 108.61955
  177.34222)

Which are significantly greater than the corresponding constraints (mapcar 'second nutrients) of

(3 70 0.8 12 5 1.8 2.7 18 75)

Unless I've miscalculated something, it appears that the solution currently produced by this library is not optimal.

Thank you so much for creating this library. I would very much like to use it. Hopefully I've just miscalculated something, or the issue is something not too difficult to fix.

Redundand constraint leads to wrong solution

;; Gives:  
;; w=1.5   x=1.0 y=0.5,
;; but y can't be 0.5  because of (>= y 1.0)
(lp:with-solved-problem
    ((min (= w (+ x y)))
     (>= x 1.0)
     (>= y 1.0)
     ;; If the first two inequalities hold, this is redundant constraint
     ;; however, it leads to wrong solution
     #-nil (>= (+ x (* 2.0 y)) 2.0))
  (format t "w=~a   x=~a y=~a~%" w x y)
  (values))

The effect is like it just takes 3rd constraint and throws away 2nd.

UNBOUNDED-PROBLEM-ERROR raised for bounded problem.

Hi,

I'm trying to solve the Stigler Diet problem using this library. For the definition of the problem as well as the code in Python, please see https://developers.google.com/optimization/lp/glop

I'm using the 2019-11-30 version of this library from Quicklisp.

When I try to run the code below, the error UNBOUNDED-PROBLEM-ERROR is signaled. Am I defining the problem incorrectly?

(defparameter data
  '(("Wheat Flour (Enriched)" "10 lb." 36 44.7 1411 2 365 0 55.4 33.3 441 0)
    ("Macaroni" "1 lb." 14.1 11.6 418 0.7 54 0 3.2 1.9 68 0)
    ("Wheat Cereal (Enriched)" "28 oz." 24.2 11.8 377 14.4 175 0 14.4 8.8 114 0)
    ("Corn Flakes" "8 oz." 7.1 11.4 252 0.1 56 0 13.5 2.3 68 0)
    ("Corn Meal" "1 lb." 4.6 36.0 897 1.7 99 30.9 17.4 7.9 106 0)
    ("Hominy Grits" "24 oz." 8.5 28.6 680 0.8 80 0 10.6 1.6 110 0)
    ("Rice" "1 lb." 7.5 21.2 460 0.6 41 0 2 4.8 60 0)
    ("Rolled Oats" "1 lb." 7.1 25.3 907 5.1 341 0 37.1 8.9 64 0)
    ("White Bread (Enriched)" "1 lb." 7.9 15.0 488 2.5 115 0 13.8 8.5 126 0)
    ("Whole Wheat Bread" "1 lb." 9.1 12.2 484 2.7 125 0 13.9 6.4 160 0)
    ("Rye Bread" "1 lb." 9.1 12.4 439 1.1 82 0 9.9 3 66 0)
    ("Pound Cake" "1 lb." 24.8 8.0 130 0.4 31 18.9 2.8 3 17 0)
    ("Soda Crackers" "1 lb." 15.1 12.5 288 0.5 50 0 0 0 0 0)
    ("Milk" "1 qt." 11 6.1 310 10.5 18 16.8 4 16 7 177)
    ("Evaporated Milk (can)" "14.5 oz." 6.7 8.4 422 15.1 9 26 3 23.5 11 60)
    ("Butter" "1 lb." 30.8 10.8 9 0.2 3 44.2 0 0.2 2 0)
    ("Oleomargarine" "1 lb." 16.1 20.6 17 0.6 6 55.8 0.2 0 0 0)
    ("Eggs" "1 doz." 32.6 2.9 238 1.0 52 18.6 2.8 6.5 1 0)
    ("Cheese (Cheddar)" "1 lb." 24.2 7.4 448 16.4 19 28.1 0.8 10.3 4 0)
    ("Cream" "1/2 pt." 14.1 3.5 49 1.7 3 16.9 0.6 2.5 0 17)
    ("Peanut Butter" "1 lb." 17.9 15.7 661 1.0 48 0 9.6 8.1 471 0)
    ("Mayonnaise" "1/2 pt." 16.7 8.6 18 0.2 8 2.7 0.4 0.5 0 0)
    ("Crisco" "1 lb." 20.3 20.1 0 0 0 0 0 0 0 0)
    ("Lard" "1 lb." 9.8 41.7 0 0 0 0.2 0 0.5 5 0)
    ("Sirloin Steak" "1 lb." 39.6 2.9 166 0.1 34 0.2 2.1 2.9 69 0)
    ("Round Steak" "1 lb." 36.4 2.2 214 0.1 32 0.4 2.5 2.4 87 0)
    ("Rib Roast" "1 lb." 29.2 3.4 213 0.1 33 0 0 2 0 0)
    ("Chuck Roast" "1 lb." 22.6 3.6 309 0.2 46 0.4 1 4 120 0)
    ("Plate" "1 lb." 14.6 8.5 404 0.2 62 0 0.9 0 0 0)
    ("Liver (Beef)" "1 lb." 26.8 2.2 333 0.2 139 169.2 6.4 50.8 316 525)
    ("Leg of Lamb" "1 lb." 27.6 3.1 245 0.1 20 0 2.8 3.9 86 0)
    ("Lamb Chops (Rib)" "1 lb." 36.6 3.3 140 0.1 15 0 1.7 2.7 54 0)
    ("Pork Chops" "1 lb." 30.7 3.5 196 0.2 30 0 17.4 2.7 60 0)
    ("Pork Loin Roast" "1 lb." 24.2 4.4 249 0.3 37 0 18.2 3.6 79 0)
    ("Bacon" "1 lb." 25.6 10.4 152 0.2 23 0 1.8 1.8 71 0)
    ("Ham smoked" "1 lb." 27.4 6.7 212 0.2 31 0 9.9 3.3 50 0)
    ("Salt Pork" "1 lb." 16 18.8 164 0.1 26 0 1.4 1.8 0 0)
    ("Roasting Chicken" "1 lb." 30.3 1.8 184 0.1 30 0.1 0.9 1.8 68 46)
    ("Veal Cutlets" "1 lb." 42.3 1.7 156 0.1 24 0 1.4 2.4 57 0)
    ("Salmon Pink (can)" "16 oz." 13 5.8 705 6.8 45 3.5 1 4.9 209 0)
    ("Apples" "1 lb." 4.4 5.8 27 0.5 36 7.3 3.6 2.7 5 544)
    ("Bananas" "1 lb." 6.1 4.9 60 0.4 30 17.4 2.5 3.5 28 498)
    ("Lemons" "1 doz." 26 1.0 21 0.5 14 0 0.5 0 4 952)
    ("Oranges" "1 doz." 30.9 2.2 40 1.1 18 11.1 3.6 1.3 10 1998)
    ("Green Beans" "1 lb." 7.1 2.4 138 3.7 80 69 4.3 5.8 37 862)
    ("Cabbage" "1 lb." 3.7 2.6 125 4.0 36 7.2 9 4.5 26 5369)
    ("Carrots" "1 bunch" 4.7 2.7 73 2.8 43 188.5 6.1 4.3 89 608)
    ("Celery" "1 stalk" 7.3 0.9 51 3.0 23 0.9 1.4 1.4 9 313)
    ("Lettuce" "1 head" 8.2 0.4 27 1.1 22 112.4 1.8 3.4 11 449)
    ("Onions" "1 lb." 3.6 5.8 166 3.8 59 16.6 4.7 5.9 21 1184)
    ("Potatoes" "15 lb." 34 14.3 336 1.8 118 6.7 29.4 7.1 198 2522)
    ("Spinach" "1 lb." 8.1 1.1 106 0 138 918.4 5.7 13.8 33 2755)
    ("Sweet Potatoes" "1 lb." 5.1 9.6 138 2.7 54 290.7 8.4 5.4 83 1912)
    ("Peaches (can)" "No. 2 1/2" 16.8 3.7 20 0.4 10 21.5 0.5 1 31 196)
    ("Pears (can)" "No. 2 1/2" 20.4 3.0 8 0.3 8 0.8 0.8 0.8 5 81)
    ("Pineapple (can)" "No. 2 1/2" 21.3 2.4 16 0.4 8 2 2.8 0.8 7 399)
    ("Asparagus (can)" "No. 2" 27.7 0.4 33 0.3 12 16.3 1.4 2.1 17 272)
    ("Green Beans (can)" "No. 2" 10 1.0 54 2 65 53.9 1.6 4.3 32 431)
    ("Pork and Beans (can)" "16 oz." 7.1 7.5 364 4 134 3.5 8.3 7.7 56 0)
    ("Corn (can)" "No. 2" 10.4 5.2 136 0.2 16 12 1.6 2.7 42 218)
    ("Peas (can)" "No. 2" 13.8 2.3 136 0.6 45 34.9 4.9 2.5 37 370)
    ("Tomatoes (can)" "No. 2" 8.6 1.3 63 0.7 38 53.2 3.4 2.5 36 1253)
    ("Tomato Soup (can)" "10 1/2 oz." 7.6 1.6 71 0.6 43 57.9 3.5 2.4 67 862)
    ("Peaches Dried" "1 lb." 15.7 8.5 87 1.7 173 86.8 1.2 4.3 55 57)
    ("Prunes Dried" "1 lb." 9 12.8 99 2.5 154 85.7 3.9 4.3 65 257)
    ("Raisins Dried" "15 oz." 9.4 13.5 104 2.5 136 4.5 6.3 1.4 24 136)
    ("Peas Dried" "1 lb." 7.9 20.0 1367 4.2 345 2.9 28.7 18.4 162 0)
    ("Lima Beans Dried" "1 lb." 8.9 17.4 1055 3.7 459 5.1 26.9 38.2 93 0)
    ("Navy Beans Dried" "1 lb." 5.9 26.9 1691 11.4 792 0 38.4 24.6 217 0)
    ("Coffee" "1 lb." 22.4 0 0 0 0 0 4 5.1 50 0)
    ("Tea" "1/4 lb." 17.4 0 0 0 0 0 0 2.3 42 0)
    ("Cocoa" "8 oz." 8.6 8.7 237 3 72 0 2 11.9 40 0)
    ("Chocolate" "8 oz." 16.2 8.0 77 1.3 39 0 0.9 3.4 14 0)
    ("Sugar" "10 lb." 51.7 34.9 0 0 0 0 0 0 0 0)
    ("Corn Syrup" "24 oz." 13.7 14.7 0 0.5 74 0 0 0 5 0)
    ("Molasses" "18 oz." 13.6 9.0 0 10.3 244 0 1.9 7.5 146 0)
    ("Strawberry Preserves" "1 lb." 20.5 6.4 11 0.4 7 0.2 0.2 0.4 3 0)))

(defparameter nutrients
  '(("Calories (1000s)" 3)
    ("Protein (grams)" 70)
    ("Calcium (grams)" 0.8)
    ("Iron (mg)" 12)
    ("Vitamin A (1000 IU)" 5)
    ("Vitamin B1 (mg)" 1.8)
    ("Vitamin B2 (mg)" 2.7)
    ("Niacin (mg)" 18)
    ("Vitamin C (mg)" 75)))

(defparameter objective-function
 `(min (= w (+ ,@(mapcar 'make-symbol (mapcar 'first data))))))

(defparameter constraints
  (let ((constraints (list)))
    (do ((i 0 (1+ i)))
        ((= i (length nutrients)) i)
      (push `(<= ,(second (nth i nutrients))
                 (+ ,@(mapcar (lambda (d) `(* ,(nth (+ i 3) d)
                                              ,(make-symbol (first d)))) data)))
            constraints))
    constraints))

(defparameter problem
  (parse-linear-problem
   objective-function
   constraints))

(defparameter solution
  (solve-problem problem))

Error while parsing arguments to MACROLET REDUCED-COST

In the README.md example, the WITH-SOLUTION-VARIABLES and WITH-SOLVED-PROBLEM macros give an error while parsing arguments to macrolet REDUCED-COST.

Using the WITH-SOLVED-PROBLEM example:

(with-solved-problem ((max (= w (+ x (* 4 y) (* 3 z))))
                      (<= (+ (* 2 x) y) 8)
                      (<= (+ y z) 7))
  (format t "Objective value solution: ~A~%" w)
  (format t "x = ~A (reduced cost: ~A)~%" x (reduced-cost solution 'x))
  (format t "y = ~A (reduced cost: ~A)~%" y (reduced-cost solution 'y))
  (format t "z = ~A (reduced cost: ~A)~%" z (reduced-cost solution 'z)))

produces the following error:

Execution of a form compiled with errors.
Form:
(REDUCED-COST SOLUTION 'X)
Compile-time error:
during macroexpansion of (REDUCED-COST SOLUTION 'X). Use BREAK-ON-SIGNALS to
intercept.

error while parsing arguments to MACROLET REDUCED-COST:
too many elements in
(SOLUTION 'X)
to satisfy lambda list
(LINEAR-PROGRAMMING/SOLVER::VAR):
exactly 1 expected, but got 2
[Condition of type SB-INT:COMPILED-PROGRAM-ERROR]

Restarts:
0: [RETRY] Retry SLIME REPL evaluation request.
1: [*ABORT] Return to SLIME's top level.
2: [ABORT] abort thread (#<THREAD "new-repl-thread" RUNNING {10005D7533}>)

Backtrace:
0: ((LAMBDA ()))
1: (SB-INT:SIMPLE-EVAL-IN-LEXENV (WITH-SOLVED-PROBLEM ((MAX #) (<= # 8) (<= # 7)) (FORMAT T "Objective value solution: A%" W) (FORMAT T "x = ~A (reduced cost: A)%" X (REDUCED-COST SOLUTION #)) (FORMA..
2: (EVAL (WITH-SOLVED-PROBLEM ((MAX #) (<= # 8) (<= # 7)) (FORMAT T "Objective value solution: A%" W) (FORMAT T "x = ~A (reduced cost: A)%" X (REDUCED-COST SOLUTION #)) (FORMAT T "y = ~A (reduced cos..
3: (SWANK::EVAL-REGION "(with-solved-problem ((max (= w (+ x (* 4 y) (* 3 z)))) ..)
4: ((LAMBDA NIL :IN SWANK-REPL::REPL-EVAL))
5: (SWANK-REPL::TRACK-PACKAGE #<CLOSURE (LAMBDA NIL :IN SWANK-REPL::REPL-EVAL) {100296D8CB}>)
6: (SWANK::CALL-WITH-RETRY-RESTART "Retry SLIME REPL evaluation request." #<CLOSURE (LAMBDA NIL :IN SWANK-REPL::REPL-EVAL) {100296D86B}>)
7: (SWANK::CALL-WITH-BUFFER-SYNTAX NIL #<CLOSURE (LAMBDA NIL :IN SWANK-REPL::REPL-EVAL) {100296D84B}>)
8: (SWANK-REPL::REPL-EVAL "(with-solved-problem ((max (= w (+ x (* 4 y) (* 3 z)))) ..)
9: (SB-INT:SIMPLE-EVAL-IN-LEXENV (SWANK-REPL:LISTENER-EVAL "(with-solved-problem ((max (= w (+ x (* 4 y) (* 3 z)))) ..)
10: (EVAL (SWANK-REPL:LISTENER-EVAL "(with-solved-problem ((max (= w (+ x (* 4 y) (* 3 z)))) ..)
11: (SWANK:EVAL-FOR-EMACS (SWANK-REPL:LISTENER-EVAL "(with-solved-problem ((max (= w (+ x (* 4 y) (* 3 z)))) ..)
12: (SWANK::PROCESS-REQUESTS NIL)
13: ((LAMBDA NIL :IN SWANK::HANDLE-REQUESTS))
14: ((LAMBDA NIL :IN SWANK::HANDLE-REQUESTS))
15: (SWANK/SBCL::CALL-WITH-BREAK-HOOK # #<CLOSURE (LAMBDA NIL :IN SWANK::HANDLE-REQUESTS) {1001B4008B}>)
16: ((FLET SWANK/BACKEND:CALL-WITH-DEBUGGER-HOOK :IN "/home/sernamar/.emacs.d/elpa/slime-2.26.1/swank/sbcl.lisp") # #<CLOSURE (LAMBDA NIL :IN SWANK::HANDLE-REQUESTS) {1..
17: (SWANK::CALL-WITH-BINDINGS ((STANDARD-INPUT . #<SWANK/GRAY::SLIME-INPUT-STREAM {1001C06123}>)) #<CLOSURE (LAMBDA NIL :IN SWANK::HANDLE-REQUESTS) {1001B400AB}>)
18: (SWANK::HANDLE-REQUESTS #<SWANK::MULTITHREADED-CONNECTION {1004CC64C3}> NIL)
19: ((FLET SB-UNIX::BODY :IN SB-THREAD::INITIAL-THREAD-FUNCTION-TRAMPOLINE))
20: ((FLET "WITHOUT-INTERRUPTS-BODY-4" :IN SB-THREAD::INITIAL-THREAD-FUNCTION-TRAMPOLINE))
21: ((FLET SB-THREAD::WITH-MUTEX-THUNK :IN SB-THREAD::INITIAL-THREAD-FUNCTION-TRAMPOLINE))
22: ((FLET "WITHOUT-INTERRUPTS-BODY-1" :IN SB-THREAD::CALL-WITH-MUTEX))
23: (SB-THREAD::CALL-WITH-MUTEX #<CLOSURE (FLET SB-THREAD::WITH-MUTEX-THUNK :IN SB-THREAD::INITIAL-THREAD-FUNCTION-TRAMPOLINE) {7FFA3F95ED6B}> #<SB-THREAD:MUTEX "thread result lock" owner: #<SB-THREAD:THR..
24: (SB-THREAD::INITIAL-THREAD-FUNCTION-TRAMPOLINE #<SB-THREAD:THREAD "new-repl-thread" RUNNING {10005D7533}> NIL #<CLOSURE (LAMBDA NIL :IN SWANK-REPL::SPAWN-REPL-THREAD) {10005D74DB}> NIL)
25: ("foreign function: call_into_lisp")
26: ("foreign function: new_thread_trampoline")

I've taken a look at the code at:

linear-programming/src/solver.lisp

Lines 92 to 111 in 2f1dc9c

	(defmacro with-solution-variables (var-list solution &body body)
	"Evaluates the body with the variables in `var-list` bound to their values in the
	solution. Additionally, the macro `reduced-cost` is locally bound that takes a
	variable name and provides it's reduced cost."
	(once-only (solution)
	(let ((body (list `(macrolet ((reduced-cost (var)
	`(solution-reduced-cost ,',solution ',var)))
	,@body))))
	(if (typep var-list 'problem)
	(let* ((problem var-list) ;alias for readability
	(vars (problem-vars problem)))
	`(let ((,(problem-objective-var problem) (solution-objective-value ,solution))
	,@(iter (for var in-vector vars)
	(for i from 0)
	(collect `(,var (solution-variable ,solution ',var)))))
	(declare (ignorable ,(problem-objective-var problem) ,@(map 'list #'identity vars)))
	,@body))
	`(let (,@(iter (for var in-sequence var-list)
	(collect `(,var (solution-variable ,solution ',var)))))
	,@body)))))

and it seems that the macrolet REDUCE-COST expects a single variable. So the following code will work as expected:

(with-solved-problem ((max (= w (+ x (* 4 y) (* 3 z))))
                      (<= (+ (* 2 x) y) 8)
                      (<= (+ y z) 7))
  (format t "Objective value solution: ~A~%" w)
  (format t "x = ~A (reduced cost: ~A)~%" x (reduced-cost x))
  (format t "y = ~A (reduced cost: ~A)~%" y (reduced-cost y))
  (format t "z = ~A (reduced cost: ~A)~%" z (reduced-cost z)))

> constraint does not behave as expected

It seems that when adding a > constraint the actual behavior is >=

E.g

(with-solved-problem ((max (= w (+ x (* 4 y) (* 3 z))))
                      (<= (+ (* 2 x) y) 8)
                      (<= (+ y z) 7)
                      (bounds (2 z 5))
                      (> x 3)
                      (integer x))
  (format t "Objective value solution: ~A~%" w)
  (format t "x = ~A (reduced cost: ~A)~%" x (reduced-cost x))
  (format t "y = ~A (reduced cost: ~A)~%" y (reduced-cost y))
  (format t "z = ~A (reduced cost: ~A)~%" z (reduced-cost z)))

Which gives

Objective value solution: 26
x = 3 (reduced cost: 7)
y = 2 (reduced cost: 0)
z = 5 (reduced cost: 0)

Add support for text input/output of problem descriptions

• Sexp format
• Human-readable, standard format writer
• MPS format
  • reader
  • writer
• Other solver formats (eg http://lpsolve.sourceforge.net/5.1/formulate.htm#Read%20the%20model%20from%20an%20ASCII%20file)