Introducción.

21. Introducción.#

El problema de la mochila es un problema clásico dentro del mundo de la optimización. La idea es que dado una serie de paquetes, con unas ciertas carcaterísticas y un lugar donde se quieren meter esos paquetes (puede ser un contenedor, camión etc), estudiar la manera óptima de colcar los paquetes en el contenedor.

En este apartado vamos a ver diferentes versiones del problema y la forma de darlo solución.

21.1. versión 1.#

Tenemos cuatro paquetes, cada paquete tiene un valor, un peso y un volumen. Se trata de maximizar el valor que se introduce teniendo el cuenta que el contenedor tiene un límite de pesos y y de volumen.

En este caso la función objetivo seria:

\[\sum valor_{i}x_{i}\]

Donde \(x_i\) es una variable dicotómica que toma el valor 1 si es seleccionado el paquete y cero en caso contrario. Observar que lo que se pretende es maximimizar, luego lo que tenemos que hacer es poner esa función con signo negativo.

Las restricciones se pondrían de forma similar. A continuación, se muestran los valores con los que vamos a trabajar

values = [1, 2, 3, 4] # Valores de los paquetes
weights = [2, 2, 3, 3] #pesos de los paquetes
volumes = [3, 3, 2, 2] # volumen de los paquetes
n = 4 # nº total de paquetes
variables = list(range(n)) # Nos va a servir despues para definir las variables binarias
max_weight = 4 # máximo peso que admite el contenedor
max_volume = 6 # Volumen máximo que admite el contenedor

penalty_strength_weight = 20 # Parámetro de lagrange para la condición de peso
penalty_strength_volume = 20 # parámetro de lagrange para la condición de volumen

A continuación mostramos una tabla de posibles valores y señalamos en negrita la opción que vamos buscando. Observar que hay soluciones que maximizan la función de costes, pero sin embargo no es una solucción factible pues no verifica alguna de las restricciones impuestas.

x_1	x_2	x_3	x_4	total value	total weight	total volume	feasible
1	0	0	0	1	2	3	True
0	1	0	0	2	2	3	True
0	0	1	0	3	3	2	True
0	0	0	1	4	3	2	True
0	0	1	1	7	6	4	False
1	1	0	0	3	4	6	True

Comenzamos implementado el modelo de tipo BQM que vamos a utilizar

from dimod import BinaryQuadraticModel

bqm = BinaryQuadraticModel('BINARY')

bqm es una instancia de BinaryQuadraticModel y además le definimos como ‘BINARY’ para indicar que las variable \(x_i\) sólo toman valores cero o uno, como ya se ha comentado antes.

variables = [bqm.add_variable(f'x_{v}', -values[v]) for v in variables]

# Veamos las variables que hemos creado,

variables

['x_0', 'x_1', 'x_2', 'x_3']

observar que a vada variable le hemos dado un ofset o coeficiente con signo negativo, porque lo que queremos es minimizar la función de coste. Veamos a continuación el modelo que hemos creado,

bqm

BinaryQuadraticModel({'x_0': -1.0, 'x_1': -2.0, 'x_2': -3.0, 'x_3': -4.0}, {}, 0.0, 'BINARY')

Ahora añadimos la restriccióm de la capacidad de la mochila. Observar que aquí para introducir una restricción de desigualdad utilizamos la propiedad bqm.add_linear_inequality_constrain t.

slacks_weight = bqm.add_linear_inequality_constraint(
    [(x, w) for x, w in zip(variables, weights)],
    constant=-max_weight, #este sería el peso máximo permitido
    lagrange_multiplier=penalty_strength_weight, # el multiplicador de lagrange definido al principio
    label='capacity'
)

Hacemos un tratamaiento similar para la restricción del volumen.

slacks_volume = bqm.add_linear_inequality_constraint(
    [(x, v) for x, v in zip(variables, volumes)],
    constant=-max_volume,
    lagrange_multiplier=penalty_strength_volume,
    label='volume_capacity'
)

Para problemas reducidos con por ejemplo menos de 20 variables, podemos resolver el problema QUBO de forma exacta enumerando todas las posibles soluciones

from dimod import ExactSolver

response = ExactSolver().sample(bqm)
samples = response.samples()
best_solution = response.first.sample
print(response.truncate(5))

  slack_capacity_0 slack_capacity_1 slack_capacity_2 ... x_3 energy num_oc.
              1                0                0 ...   1   -4.0       1
              0                0                1 ...   1   -4.0       1
              0                0                0 ...   0   -3.0       1
              1                0                0 ...   0   -3.0       1
              0                0                1 ...   0   -3.0       1
['BINARY', 5 rows, 5 samples, 10 variables]

best_solution.items()

dict_items([('slack_capacity_0', 1), ('slack_capacity_1', 0), ('slack_capacity_2', 0), ('slack_volume_capacity_0', 1), ('slack_volume_capacity_1', 0), ('slack_volume_capacity_2', 1), ('x_0', 0), ('x_1', 0), ('x_2', 0), ('x_3', 1)])

Vamos a examinar la solución a continuación.

Las variables principales especifican el montaje de la mochila.
El primer conjunto de variables de holgura indica que el peso total de los artículos seleccionados es 1 unidad inferior a la capacidad de peso de la mochila.
El segundo conjunto de variables de holgura indica que el volumen total de los artículos seleccionados es 1 unidad inferior a la capacidad de volumen de la mochila.

print({k: v for k, v in best_solution.items() if k in variables})
print({k: v for k, v in best_solution.items() if k in {v for v, _ in slacks_weight}})
print({k: v for k, v in best_solution.items() if k in {v for v, _ in slacks_volume}})

{'x_0': 0, 'x_1': 0, 'x_2': 0, 'x_3': 1}
{'slack_capacity_0': 1, 'slack_capacity_1': 0, 'slack_capacity_2': 0}
{'slack_volume_capacity_0': 1, 'slack_volume_capacity_1': 0, 'slack_volume_capacity_2': 1}

Definimos a continuación una función que nos permita saber si se cumplen o no las restricciones del problema

def violations(sample, slack_variables, array, capacity):
    slack_total_weight = 0
    weight = 0
    for i, (v, c) in enumerate(slack_variables):
        slack_total_weight += sample[v] * c
    for i, v in enumerate(variables):
        weight += array[i] * sample[v]
    print('Total of items: ', weight)
    print('Total slack value: ', slack_total_weight)
    print('Left hand side: ', weight + slack_total_weight)
    print('Satisfied? ', weight + slack_total_weight <= capacity)

print('Weight capacity constraint:')
violations(samples[0], slacks_weight, weights, max_weight)
print('Volume capacity constraint:')
violations(samples[0], slacks_volume, volumes, max_volume)

Weight capacity constraint:
Total of items:  3
Total slack value:  1
Left hand side:  4
Satisfied?  True
Volume capacity constraint:
Total of items:  2
Total slack value:  4
Left hand side:  6
Satisfied?  True

21.2. Version 2 (con CQM).#

El problema es similar al mostrado anteriormente, pero en esta ocasión utilizo una instancia de la clase CQM. Este tipo de solución también lo podemos ver en https://github.com/dwave-examples/knapsack .

from dimod import Binary, CQM, quicksum
from dwave.system import LeapHybridCQMSampler
import itertools

values = [34, 25,78,21,64]
weights = [3,5,9,4,2]
W = 10
n= len(values)

 #Creamos las variable binarias
x = [Binary(i)  for i in range(n)]

 #Creamos las variable binarias
x = [Binary(i)  for i in range(n)]

 cqm = CQM()
# Añadimos la función objetivo
cqm.set_objective(quicksum(-values[i]*x[i] for i in range(n)))

# Añadimos las restricciones
cqm.add_constraint(quicksum(weights[i]*x[i] for i in range(n)) <= W, label ='peso maximo')
cqm.add_constraint(quicksum(x[i] for i in range(n))<=2, label ='items máximo ')

'items máximo '

# submitimos el sampler
sampler =  LeapHybridCQMSampler()
sampleset = sampler.sample_cqm(cqm)

21.3. Versión 3 (Con otras restricciones).#

En este problema se trata de empaquetar un camión de reparto. En este ejemplo, consideramos un escenario simple en el que tenemos un solo camión y un solo tamaño de caja. Queremos elegir paquetes para cargar en el camión de modo que se maximice la cantidad de paquetes con estado de envío prioritario y se minimice el tiempo de tránsito total. No se debe exceder la capacidad de peso y tamaño del camión.

from dimod import ConstrainedQuadraticModel, Binaries, quicksum 
from dwave.system import LeapHybridCQMSampler 
import random 
import numpy as np

# Problem set up
num_packages = 300

# Prioridad de cada paquete, 3 = high priority, 1 = low priority
priority = random.choices((1, 2, 3), k=num_packages)

# Número de días transcurridos desde que se encargó cada paquete (a más días, mayor prioridad)
days_since_order = random.choices((0, 1, 2, 3), k=num_packages)

# Weight of each package
weight = random.choices(range(1, 101), k=num_packages)

# Peso máximo y número de bultos que puede transportar el camión
max_weight = 3000  
max_packages = 100

# Pesos para cada sumando de la función objetivo
obj_weight_priority = 1.0 
obj_weight_days = 1.0

# Build the CQM
cqm = ConstrainedQuadraticModel()

# Create the binary variables
bin_variables = list(Binaries(range(num_packages)))

# Veamos cómo es una variable
bin_variables[0]

BinaryQuadraticModel({0: 1.0}, {}, 0.0, 'BINARY')

# quicksum() se utiliza para indicar que es una suma

# ----------------- Objective functions ----------------- 

# 1. Maximize priority shipping packages
objective1 = -obj_weight_priority * quicksum(priority[i] * bin_variables[i] 
                for i in range(num_packages))

# 2. Minimize customers wait time
objective2 = -obj_weight_days * quicksum(days_since_order[i] * bin_variables[i]
                for i in range(num_packages))

# Add the objectives to the CQM
cqm.set_objective(objective1 + objective2)

# ----------------- Constraints ----------------- 

# Add the maximum capacity constraint
cqm.add_constraint(quicksum(weight[i] * bin_variables[i] for i in
range(num_packages)) <= max_weight, label='max_weight')

# Add the maximum package (or truck size) constraint
cqm.add_constraint(quicksum(bin_variables[i] for i in range(num_packages)) 
                    == max_packages, label='max_packages')

'max_packages'

# -----------------  Submit to the CQM sampler -----------------
cqm_sampler = LeapHybridCQMSampler() 
sampleset = cqm_sampler.sample_cqm(cqm, label='Truck Packing Demo')

feasible_sampleset = sampleset.filter(lambda d: d.is_feasible)

if not len(feasible_sampleset): 
    print("\nNo feasible solution found.\n")

else: 
    first_feasible_sol = feasible_sampleset.first.sample
    print(first_feasible_sol)

{0: 0.0, 1: 0.0, 2: 0.0, 3: 0.0, 4: 0.0, 5: 0.0, 6: 0.0, 7: 1.0, 8: 0.0, 9: 1.0, 10: 0.0, 11: 0.0, 12: 0.0, 13: 0.0, 14: 0.0, 15: 0.0, 16: 0.0, 17: 1.0, 18: 1.0, 19: 0.0, 20: 0.0, 21: 0.0, 22: 0.0, 23: 0.0, 24: 1.0, 25: 1.0, 26: 1.0, 27: 1.0, 28: 1.0, 29: 0.0, 30: 0.0, 31: 0.0, 32: 0.0, 33: 0.0, 34: 1.0, 35: 0.0, 36: 1.0, 37: 0.0, 38: 0.0, 39: 1.0, 40: 0.0, 41: 1.0, 42: 1.0, 43: 0.0, 44: 1.0, 45: 1.0, 46: 1.0, 47: 0.0, 48: 0.0, 49: 0.0, 50: 1.0, 51: 0.0, 52: 0.0, 53: 1.0, 54: 1.0, 55: 0.0, 56: 0.0, 57: 0.0, 58: 0.0, 59: 1.0, 60: 0.0, 61: 0.0, 62: 1.0, 63: 0.0, 64: 0.0, 65: 0.0, 66: 1.0, 67: 1.0, 68: 1.0, 69: 0.0, 70: 0.0, 71: 0.0, 72: 0.0, 73: 0.0, 74: 1.0, 75: 1.0, 76: 0.0, 77: 0.0, 78: 0.0, 79: 1.0, 80: 1.0, 81: 0.0, 82: 0.0, 83: 0.0, 84: 0.0, 85: 0.0, 86: 0.0, 87: 1.0, 88: 1.0, 89: 0.0, 90: 0.0, 91: 1.0, 92: 1.0, 93: 0.0, 94: 0.0, 95: 0.0, 96: 1.0, 97: 1.0, 98: 0.0, 99: 0.0, 100: 0.0, 101: 0.0, 102: 0.0, 103: 0.0, 104: 0.0, 105: 1.0, 106: 0.0, 107: 0.0, 108: 0.0, 109: 0.0, 110: 0.0, 111: 0.0, 112: 0.0, 113: 0.0, 114: 0.0, 115: 0.0, 116: 0.0, 117: 0.0, 118: 0.0, 119: 1.0, 120: 0.0, 121: 1.0, 122: 0.0, 123: 0.0, 124: 0.0, 125: 0.0, 126: 0.0, 127: 1.0, 128: 0.0, 129: 0.0, 130: 0.0, 131: 0.0, 132: 1.0, 133: 0.0, 134: 1.0, 135: 0.0, 136: 0.0, 137: 0.0, 138: 0.0, 139: 0.0, 140: 1.0, 141: 0.0, 142: 0.0, 143: 1.0, 144: 1.0, 145: 0.0, 146: 0.0, 147: 1.0, 148: 0.0, 149: 0.0, 150: 0.0, 151: 1.0, 152: 0.0, 153: 1.0, 154: 0.0, 155: 0.0, 156: 0.0, 157: 0.0, 158: 0.0, 159: 0.0, 160: 0.0, 161: 0.0, 162: 0.0, 163: 1.0, 164: 1.0, 165: 0.0, 166: 0.0, 167: 0.0, 168: 0.0, 169: 0.0, 170: 0.0, 171: 1.0, 172: 0.0, 173: 0.0, 174: 0.0, 175: 1.0, 176: 0.0, 177: 1.0, 178: 0.0, 179: 0.0, 180: 1.0, 181: 0.0, 182: 0.0, 183: 0.0, 184: 1.0, 185: 0.0, 186: 1.0, 187: 0.0, 188: 0.0, 189: 0.0, 190: 1.0, 191: 0.0, 192: 0.0, 193: 0.0, 194: 0.0, 195: 0.0, 196: 0.0, 197: 0.0, 198: 1.0, 199: 0.0, 200: 0.0, 201: 1.0, 202: 1.0, 203: 0.0, 204: 0.0, 205: 0.0, 206: 0.0, 207: 0.0, 208: 1.0, 209: 0.0, 210: 1.0, 211: 0.0, 212: 0.0, 213: 1.0, 214: 0.0, 215: 0.0, 216: 0.0, 217: 1.0, 218: 1.0, 219: 0.0, 220: 1.0, 221: 0.0, 222: 0.0, 223: 0.0, 224: 1.0, 225: 0.0, 226: 0.0, 227: 1.0, 228: 0.0, 229: 1.0, 230: 0.0, 231: 1.0, 232: 0.0, 233: 1.0, 234: 0.0, 235: 0.0, 236: 0.0, 237: 1.0, 238: 0.0, 239: 0.0, 240: 0.0, 241: 0.0, 242: 1.0, 243: 1.0, 244: 0.0, 245: 0.0, 246: 1.0, 247: 0.0, 248: 0.0, 249: 0.0, 250: 1.0, 251: 0.0, 252: 1.0, 253: 1.0, 254: 0.0, 255: 0.0, 256: 0.0, 257: 0.0, 258: 1.0, 259: 1.0, 260: 1.0, 261: 1.0, 262: 1.0, 263: 1.0, 264: 0.0, 265: 0.0, 266: 0.0, 267: 0.0, 268: 1.0, 269: 1.0, 270: 1.0, 271: 0.0, 272: 0.0, 273: 0.0, 274: 0.0, 275: 1.0, 276: 1.0, 277: 0.0, 278: 1.0, 279: 0.0, 280: 0.0, 281: 1.0, 282: 1.0, 283: 1.0, 284: 1.0, 285: 1.0, 286: 1.0, 287: 1.0, 288: 1.0, 289: 0.0, 290: 0.0, 291: 1.0, 292: 0.0, 293: 0.0, 294: 0.0, 295: 1.0, 296: 0.0, 297: 1.0, 298: 0.0, 299: 0.0}

# ----------------- Process the results -----------------
feasible_sampleset = sampleset.filter(lambda d: d.is_feasible)

if not len(feasible_sampleset): 
    print("\nNo feasible solution found.\n")

else: 
    first_feasible_sol = feasible_sampleset.first.sample

    # Characterize the problem
    problem_array = np.zeros((3, 4)).astype(int) 
    for i in range(num_packages):
        problem_array[-1 * (priority[i]-3)][-1 * (days_since_order[i] - 3)] += 1

    print("\n****************** PROBLEM ******************\n") 
    print('Days since order was placed')
    print('{:>5s}{:>5s}{:>5s}{:>5s}{:>5s}'.format('Priority |',
             '3', '2', '1', '0')) 
    print('-' * 40)

    for i in range(3): 
        print('{:>5s}{:>10s}{:>5s}{:>5s}{:>5s}'.format(str(-1*(i- 3)), 
                str(problem_array[i][0]), str(problem_array[i][1]),
                str(problem_array[i][2]), str(problem_array[i][3])))

    # Calculate number of packages with each priority and number of days since
    # order in the solution
    chosen = [i for i in first_feasible_sol if first_feasible_sol[i] == 1]
    total_weight = quicksum(weight[i] for i in chosen) 
    solution_priorities = [priority[i] for i in chosen] 
    solution_days_since_order = [days_since_order[i] for i in chosen]

    # Characterize the solution 
    # Packages with higher priority (upper row) and the most days since the
    # order (left most column) should be prioritized in the selection
    results_array = np.zeros((3, 4)).astype(int) 
    for i in chosen:
        results_array[-1 * (priority[i]-3)][-1 * (days_since_order[i] - 3)] += 1

    print("\n****************** SOLUTION ******************\n") 
    print('Days since order was placed')
    print('{:>5s}{:>5s}{:>5s}{:>5s}{:>5s}'.format('Priority |',
            '3', '2', '1', '0')) 
    print('-' * 40)

    for i in range(3): 
        print('{:>5s}{:>10s}{:>5s}{:>5s}{:>5s}'.format(str(-1*(i - 3)),
                str(results_array[i][0]), str(results_array[i][1]),
                str(results_array[i][2]), str(results_array[i][3])))

    print(("\nTotal number of selected items: {}".format(len(chosen))))
    print("Total weight of selected items: {}".format(total_weight))

****************** PROBLEM ******************

Days since order was placed
Priority |    3    2    1    0
----------------------------------------
    3        21   24   23   16
    2        30   22   27   31
    1        19   35   25   27

****************** SOLUTION ******************

Days since order was placed
Priority |    3    2    1    0
----------------------------------------
    3        20   15    7    3
    2        23   12    3    0
    1        12    5    0    0

Total number of selected items: 100
Total weight of selected items: 2998

Introducción.

Contents

21. Introducción.#

21.1. versión 1.#

21.2. Version 2 (con CQM).#

21.3. Versión 3 (Con otras restricciones).#