<jupyter_start><jupyter_text>Knapsack Cover Inequalities Consider the 0-1 knapsack Set $K := \Big\{\{0, 1\}^{n} : \sum_{j = 1}^n w_j x_j \le b \Big\}, w_j > 0, \forall j \in N = \{1, 2, ..., n\}, b > 0$**Definition (Cover)** A set $C \subseteq N$ is a *cover* if $\sum_{j \in C} w_j > b$. A cover is minimal if $\sum_{j \in C \backslash \{j\}} w_j \le b, \forall j \in C$. **Prop** If $C \subseteq N$ is a cover for $K$, then $\sum_{j \in C} \le \vert C \vert - 1$ is valid for $K$. The alternate reformulation of the Knapsack set is in terms of cover inequalities. In this tutorial, we will see how to generate a minimal cover, how to lift it to make it facet defining, and solve the resulting Knapsack problem.<jupyter_code># First we import some packages
import numpy as np
from gurobipy import *
import time<jupyter_output><empty_output><jupyter_text>Generating problem instance<jupyter_code>def generateData(n): # n is the number of items
    np.random.seed(30) # Setting a seed value
    a = np.random.randint(1000, size=(n)) # Generating value of items 
    w = np.random.randint(1000, size=(n))# Generating weight of items 
    b = np.sum(w)*0.40 # Total weight of Knapsack
    return a, w, b<jupyter_output><empty_output><jupyter_text>Simply solving the Knapsack problem using Gurobi<jupyter_code>def solveKnapsackUsingGurobi(n, a, w, b, verbose = 1): # no. of items, value of items, weight of items, and Knapsack weight
    m = Model() # Preparing the model
    x = {i : m.addVar(lb=0, vtype=GRB.BINARY) for i in range(n)} # whether the item should be included in the Knapsack
    m.addConstr(sum([w[j] * x[j] for j in range(n)]) <= b) # Total weight should be less than the Knapsack capacity
    m.setObjective(sum([a[j] * x[j] for j in range(n)]), sense = GRB.MAXIMIZE) # Maximize the total value of the Knapsack
    m.Params.OutputFlag = 0; m.optimize()
    if m.status == 2:  
        if verbose == 1:            
            print('*********************** Gurobi ***************************')
            print('Items selected: ', {k for k in x if x[k].x > 0.4})
            print('Total value: ', m.objVal)
            print('**********************************************************')
        return {k for k in x if x[k].x > 0.4}, m.objVal
    else:
        print('Infeasible model')
        return {}, -float("inf")<jupyter_output><empty_output><jupyter_text>Let's start preparing the initial model<jupyter_code>def preparePartialModel(n, a, w, b):
    m = Model()
    x = {i : m.addVar(lb=0.0, ub = 1.0, vtype=GRB.CONTINUOUS, name = str(i)) for i in range(n)}
    m.addConstr(sum([w[j] * x[j] for j in range(n)]) <= b)
    c, zeta = generateCover(n, a, w, b, {k:0 for k in range(n)})      
        
    m.addConstr(sum([x[j] for j in c]) <= len(c)-1)
    m.setObjective(sum([a[j] * x[j] for j in range(n)]), sense = GRB.MAXIMIZE)
    m.Params.OutputFlag = 0; m.Params.IntegralityFocus = 1
    return m<jupyter_output><empty_output><jupyter_text>Generating a cover/Seperation problemTo use cover inequalities in a cutting plane scheme, one is faced with the separation problem, that is, given a vector $\bar{x} \in \left[0, 1\right]^n$, find a cover inequality for $K$ that is violated by $\bar{x}$, or show that none exists. For this purpose, we solve the following problem:$\zeta = \min\{\sum_{j \in C} (1-\bar{x}_j): C \text{ is a cover for } K\}$- If $\zeta \ge 1$, then $\bar{x}$ satisfies all the cover inequalities for $K$. - If $\zeta < 1$, then an optimal over to the above problem yields the violated minimal cover.One can also write the above prblem as the following binary program. Let $z_j = 1$, if $j$ is in the violated cover $C$.$$	\begin{align}	& \zeta = \underset{z}{\text{minimize}}	& & \sum_{j=1}^n (1-\bar{x}_j) z_j \\	& \text{subject to}	& & \sum_{j= 1}^n w_j z_j \ge b + 1\\	& & & z \in \{0,1\}^n	\end{align}$$<jupyter_code>def generateCover(n, a, w, b, xbar): # no. of items, value of items, weight of items, Knapsack weight, and xbar
     m1 = Model()
     z = {i:m1.addVar(lb=0.0, vtype=GRB.BINARY) for i in range(n)}
     m1.addConstr(sum([w[j] * z[j] for j in range(n)]) >= b + 1)
     m1.setObjective(sum([(1-xbar[j]) * z[j] for j in range(n)]), sense = GRB.MINIMIZE)
     m1.Params.OutputFlag = 0; m1.optimize()
     return [j for j in range(n) if z[j].x > 0.4], m1.objVal # return the cover and the value of zeta<jupyter_output><empty_output><jupyter_text>Lifting This is esentially a procedure to convert a minimal cover inequality $\sum_{j \in C} x_j \le \vert C \vert - 1$ into a facet defining inequality for $K$. Sequential Lifting- Choose an ordering $j_1, j_2, ..., j_l$ of the indices in $N \backslash C$. Set t = 1.- The valid inequality $\sum_{i = 1}^{t-1} \alpha_{j_i} x_{j_i} + \sum_{j \in C} x_j \le \vert C \vert - 1$  has been obtained so far.- To calculate the largest value of $\alpha_{j_t}$ for which the inequality$\alpha_{j_t} x_t  + \sum_{i = 1}^{t-1} \alpha_{j_i} x_{j_i} + \sum_{j \in C} x_j \le \vert C \vert -1$ is valid, solve the following Knapsack problem:$$	\begin{align}	& \xi_t = \underset{x}{\text{maximize}}	& & \sum_{i = 1}^{t-1} \alpha_{j_i} x_{j_i} + \sum_{j \in C} x_j \\	& \text{subject to}	& & \sum_{i = 1}^{t-1} \alpha_{j_i} x_{j_i} + \sum_{j \in C} w_j x_j \le b - w_{j_t}\\	& & & x \in \{0,1\}^{\vert C \vert + t - 1}	\end{align}$$- Set  $\alpha_{j_t} = \vert C \vert - 1 - \xi_t$.- Stop if $t = r$.We remark that different orderings of $N \backslash C$ may produce different lifted inequalities. Furthermore, not all possible liftings can be derived from the above procedure.<jupyter_code>def sequentialLifting(c, n, a, w, b): # Cover set, no. of items, value of items, weight of items, and Knapsack weight.
    alpha = {k:0.0 for k in range(n) if k not in c} # coefficients of x_j, j not in cover
    m2 = Model()
    x_dash = {i:m2.addVar(lb=0.0, vtype=GRB.BINARY) for i in c}
    toBechecked = list(alpha.keys()); checked = [] # choosing the ordering of indices one-by-one and finding its coeffcient
    if len(toBechecked) != 0:        
        i = toBechecked.pop(0); 
        x_dash[i] = m2.addVar(lb=0.0, vtype=GRB.BINARY) # including a variable.
        m2.update()
        tmpc = sum([x_dash[j] * w[j] for j in c]) + sum([alpha[j]*x_dash[j] for j in checked]) # constraint of already in the Knapsack
        tmpo = sum([x_dash[j]  for j in c]) + sum([alpha[j]*x_dash[j] for j in checked])
        m2.setObjective(tmpo, sense = GRB.MAXIMIZE)
        m2.addConstr(tmpc <= b - w[i])
        m2.Params.OutputFlag = 0; m2.optimize()       
        if len(c) - 1 - m2.objVal >= 1: # Checking if the current item can be included in the cover.
            alpha[i] = len(c) - 1 - m2.objVal
            checked.append(i)
    return alpha # return coefficients of x_j in the lifted inequality.<jupyter_output><empty_output><jupyter_text>Balas' lifiting of minimal cover inequalitiesLet $C$ be a minimal cover for $K$, and let $\sum_{j \in C} x_j + \sum_{j \in N \backslash C} w_j x_j \le \vert C \vert - 1$ be a lifting of the cover inequality associated with $C$. Up to permuting the indices, assume that C = \{1, ..., t\} and $a_1 \ge a_2 \ge ... \ge a_t$. Let $\mu_0 =  0$, and $\mu_h = \sum_{l = 0}^h w_j,$ for $h = 1, ..., t$	. Let $\lambda = \mu_t - b$ (note that $\lambda > 0$).If above inequality defines a facet of $conv(K)$, then the following hold for every $j \in N \backslash C$.- If, for some $h$, $\mu_h ≤ w_j ≤ \mu_{h+1} − \lambda$, then $\alpha_j = h$.- If, for some $h$, $\mu_{h+1} − \lambda < w_j < \mu_{h+1}$, then $h ≤ \alpha_j ≤ h + 1$ (i.e., select $\alpha_j = h+1$).<jupyter_code>def BalasLifting(c, n, a, w, b): # Cover set, no. of items, value of items, weight of items, and Knapsack weight.
    orderCoverWeightsList = sorted(c, key=lambda item: w[item], reverse =True)
    alpha = {k:0.0 for k in range(n) if k not in c} # coefficients of x_j, j not in cover
    mu = {0:0.0}   
    ind = 0; cum = 0
    for k in orderCoverWeightsList:
        if ind == 0:
            mu[ind+1] = w[k]
            ind += 1; cum += w[k]
        else:
            mu[ind+1] = w[k] + cum
            ind += 1; cum += w[k]
    lamb = mu[len(mu)-1] - b
    
    
    for j in alpha:
        h = [h for h in mu if h < len(mu)-1 and mu[h] <= w[j] <= mu[h+1] - lamb]
        if len(h) != 0:
            alpha[j] = h[0]
        else:
            h = [h for h in mu if h < len(mu)-1 and mu[h] -lamb <= w[j] <= mu[h]]
            alpha[j] = h[-1]
    return alpha # return coefficients of x_j in the lifted inequality.<jupyter_output><empty_output><jupyter_text>Solving Knapsack problem using cover inequalitiesPutting above peices together.<jupyter_code>def solveKnapsackUsingCoverInequalities(n, a, w, b, lifting = 'sequential', verbose = 1):
    m = preparePartialModel(n, a, w, b)
    m.optimize()
    if m.status == 2:
        x = {i:m.getVarByName(str(i)).x for i in range(n)}   
        c, zeta = generateCover(n, a, w, b, x)        
    else:    
        print('Infeasible model')       
    alpha = []
   
    while zeta < 1:
        m.addConstr(sum([m.getVarByName(str(i)) for i in c]) + sum([m.getVarByName(str(i)) * alpha[i] for i in alpha]) <= len(c) - 1)
        m.optimize()
        x = {i:m.getVarByName(str(i)).x for i in range(n)}   
        c, zeta = generateCover(n, a, w, b, x) 
        #print('Cover', sum([w[k] for k in c])-min([w[k] for k in c]) <= b)
        if lifting == 'sequential':
            alpha = sequentialLifting(c, n, a, w, b)
        elif lifting == 'Balas':
            alpha = BalasLifting(c, n, a, w, b)
        elif lifting == 'None':
            alpha = []
        else:
            print('lifting is unknown.. terminating')
            return 'NA'
        zeta = round(zeta, 4)
        #print(c, zeta, alpha)
        
    if verbose == 1:
        print('******************* Cover Inequality *********************')
        print('Items selected: ', {i:m.getVarByName(str(i)).x for i in range(n) if m.getVarByName(str(i)).x > 0.4} )
        print('Total value: ', sum([a[k] for k in range(n) if round(m.getVarByName(str(k)).x) > 0.4]))
        print('**********************************************************')
    return {i:m.getVarByName(str(i)).x for i in range(n) if m.getVarByName(str(i)).x > 0.4}, sum([a[k] for k in range(n) if round(m.getVarByName(str(k)).x) > 0.4]), len(m.getConstrs()) - 1<jupyter_output><empty_output><jupyter_text>Solve the models and print results<jupyter_code>import pandas as pd
n = 1000 # No. of items
results = []; results.append(['Total items', n])
a, w, b = generateData(n)
start = time.time()
xg, objg = solveKnapsackUsingGurobi(n, a, w, b, verbose = 0)
#print('Gurobi model took', time.time()-start); start = time.time()
timeF = time.time()-start
results.append(['No. of items selected', len(xg)])
results.append(['Value of items selected', objg])
results.append(['No. of covers generated', float("NAN")])
results.append(['Computational time', timeF])

for lift in ['None', 'sequential', 'Balas']:    
    start = time.time()
    xc, objc, nc = solveKnapsackUsingCoverInequalities(n, a, w, b, lifting = lift, verbose = 0)
    timeF = time.time()-start
    #print('CG model took', time.time()- start)
    results[0].append(n)
    results[1].append(len(xc))
    results[2].append(objc)
    results[3].append(nc)
    results[4].append(timeF)

# Create the pandas DataFrame
pd.options.display.float_format = "{:,.2f}".format
df = pd.DataFrame(results, columns = ['Particulars', 'IP with Gurobi', 'Knapsack cover with no lifting', 'Knapsack cover with sequential lifting', 'Knapsack cover with Balas lifting'])
display(df)<jupyter_output>Academic license - for non-commercial use only - expires 2021-09-17
Using license file C:\Users\<NAME>\gurobi.lic
{'fraction_non_alphanumeric': 0.1729932611106372, 'fraction_numerical': 0.01484261707753988, 'mean_word_length': 3.0764951321279557, 'pattern_counts': {'":': 0, '<': 53, '<?xml version=': 0, '>': 50, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '9', 'dir': 'jupyter-structured-clean-dedup', 'id': '136510', 'license': 'MIT', 'n_tokens_mistral': 4394, 'n_tokens_neox': 4231, 'n_words': 1402, 'path': 'Knapsack Cover Inequalities/Knapsack Cover Inequalities notebook.ipynb', 'repo_name': 'prameshk/Discrete-Optimization', 'size': '11730'}