Ticket #6749: knapsack.2.patch

sage/numerical/knapsack.py

@@ -6 +6 @@
 following knapsack problems are implemented:
 
 - Solving the subset sum problem for super-increasing sequences.
-
+- General case using Linear Programming     
 
 AUTHORS:
 
 - Minh Van Nguyen (2009-04): initial version
+- Nathann Cohen (2009-08) : Linear Programming version         
 
+Definition of Knapsack problems
+---------------------------------
 
-EXAMPLES:
+You have already had a knapsack problem, so you should know, but in case you
+do not, a knapsack problem is what happens when you have hundred of items to 
+put into a bag which is too small for all of them. 
+
+When you formally write it, here is your problem :
+* Your bag can contain a weight of at most W
+* Each item `i` you have has a weight `w_i`
+* Each item `i` has a usefulness of `u_i`
+
+You then want to maximize the usefulness of the items you will store into 
+your bag, while keeping sure the weight of the bag will not go over W
+    
+As a linear program, this problem can be represented this way 
+( if you define `b_i` as the binary variable indicating whether 
+the item `i` is to be included in your bag ) :
+
+.. MATH::
+    \mbox{Maximize : }\sum_ib_iu_i\\
+    \mbox{Such that : }
+    \sum_ib_iw_i\leq W\\
+    \forall i, b_i \mbox{ binary variable}
+    
+( for more informations, cf. http://en.wikipedia.org/wiki/Knapsack_problem )
+
+Examples
+--------------------
+    
+If your knapsack problem is composed of three items (weight, value) 
+defined by (1,2), (1.5,1), (0.5,3), and a bag of maximum weight 2,
+you can easily solve it this way ::
+
+    sage: from sage.numerical.knapsack import knapsack
+    sage: knapsack( [(1,2), (1.5,1), (0.5,3)], max=2) # optional - requires Glpk or COIN-OR/CBC
+    [5.0, [(1, 2), (0.500000000000000, 3)]]
+
+Super-increasing sequences :
 
 We can test for whether or not a sequence is super-increasing::
 
@@ -33 +71 @@
     sage: L = [1, 2, 5, 21, 69, 189, 376, 919]
     sage: Superincreasing(L).subset_sum(98)
     [69, 21, 5, 2, 1]
+
+
+Functions and methods
+---------------------
 """
 
 #*****************************************************************************
@@ -509 +551 @@
             return candidates
         else:
             return []
+
+def knapsack(seq,binary=True,max=1, value_only=False):
+    r"""
+    Solves the knapsack problem
+
+
+    Knapsack problems :
+
+    You have already had a knapsack problem, so you should know, 
+    but in case you do not, a knapsack problem is what happens 
+    when you have hundred of items to put into a bag which is 
+    too small for all of them. 
+
+    When you formally write it, here is your problem :
+    * Your bag can contain a weight of at most W
+    * Each item `i` you have has a weight `w_i`
+    * Each item `i` has a usefulness of `u_i`
+
+    You then want to maximize the usefulness of the items you 
+    will store into your bag, while keeping sure the weight of 
+    the bag will not go over W.
+    
+    As a linear program, this problem can be represented this way 
+    ( if you define `b_i` as the binary variable indicating whether 
+    the item `i` is to be included in your bag ) :
+
+    .. MATH::
+        \mbox{Maximize : }\sum_ib_iu_i\\
+        \mbox{Such that : } 
+        \sum_ib_iw_i\leq W\\
+        \forall i, b_i \mbox{ binary variable}\\
+    
+    ( for more informations, 
+    cf. http://en.wikipedia.org/wiki/Knapsack_problem )
+
+    EXAMPLE :
+    
+    If your knapsack problem is composed of three 
+    items (weight, value) defined by (1,2), (1.5,1), (0.5,3), 
+    and a bag of maximum weight 2, you can easily solve it this way ::
+
+        sage: from sage.numerical.knapsack import knapsack
+        sage: knapsack( [(1,2), (1.5,1), (0.5,3)], max=2) # optional - requires Glpk or COIN-OR/CBC
+        [5.0, [(1, 2), (0.500000000000000, 3)]]
+
+        sage: knapsack( [(1,2), (1.5,1), (0.5,3)], max=2, value_only=True) # optional - requires Glpk or COIN-OR/CBC
+        5.0
+
+    In the case where all the values (usefulness) of the items 
+    are equal to one, you do not need embarass yourself with 
+    the second values, and you can just type for items 
+    `(1,1), (1.5,1), (0.5,1)` the command ::
+
+        sage: from sage.numerical.knapsack import knapsack
+        sage: knapsack([1,1.5,0.5], max=2, value_only=True) # optional - requires Glpk or COIN-OR/CBC
+        2.0
+
+    INPUT:
+    
+    - ``seq`` :  Two different possible types :
+                  - A sequence of pairs (weight, value)
+                  - A sequence of reals (a value of 1 is assumed)
+    
+    - ``binary`` : When set to True, an item can be taken 0 or 1 time.
+                   When set to False, an item can be taken any amount of
+                   times ( while staying integer and positive )
+
+    - ``max`` : Maximum admissible weight
+
+    - ``value_only`` : When set to True, only the maximum useful 
+                         value is returned
+                       When set to False, both the maximum useful 
+                         value and an assignment are returned
+
+    OUTPUT:
+    
+    If value_only is set to True, only the maximum useful value 
+    is returned.
+    Else ( default ), the function returns a pair ``[value,list]``, 
+    where ``list`` can be of two types according to the type of ``seq`` :
+
+        - A list of pairs (w_i, u_i), for each object i occurring 
+            in the solution.
+        - A list of reals where each real is repeated the number 
+            of times it is taken into the solution.
+    """
+    reals=not isinstance(seq[0],tuple)
+    if reals:
+        seq=[(x,1) for x in seq]
+
+    from sage.numerical.mip import MIP
+    p=MIP(sense=1)
+    p.setobj(dict([(i,seq[i][1]) for i in range(len(seq))]))
+    p.addconstraint(dict([(i,seq[i][0]) for i in range(len(seq))]),max=max)
+
+    if binary:
+        for i in range(len(seq)):
+            p.setbinary(i)
+    else:
+        for i in range(len(seq)):
+            p.setinteger(i)
+    
+    if value_only:
+        return p.solve(objective_only=True)
+    
+    else:
+        tmp=p.solve()
+
+        if reals:
+            val=[x[0] for c in [ [seq[a]]*int(round(b)) for (a,b) in tmp[1].items() if b>0] for x in c]
+        else:
+            val=[x for c in [ [seq[a]]*int(round(b)) for (a,b) in tmp[1].items() if b>0] for x in c]
+
+        return [tmp[0],val]