Programación lineal mixta con Google OR-Tools

En artículos anteriores hemos mencionado la diferencia existente entre programación lineal (PL) y programación lineal entera (PLE). Recordamos entonces que, cuando un modelo presenta todas sus variables enteras, se denomina puro. En caso contrario, cuando utiliza una combinación de variables enteras y continuas, se denomina mixto, constituyendo un modelo de programación lineal mixta.

En materia de optimización lineal, la programación lineal mixta, lógicamente, aborda la mayor cantidad de casos de aplicación práctica. En la medida en que se consideren la mayor cantidad de variables que representen mediante un modelo de optimización, la realidad; es potencialmente adecuado, hacer uso de modelos de programación mixta.

Por ejemplo, pensemos en la producción de televisores, es posible que en un modelo lineal queramos representar por medio de una variable entera, la cantidad de unidades producidas; ahora bien, también es posible que queramos representar por medio de una variable continua, el tiempo empleado en la fabricación. En este escenario, requerimos de programación lineal mixta.

El objetivo de este artículo consiste en utilizar las librerías del software Google OR-Tools para abordar problemas de programación lineal mixta (optimización lineal mixta). 

OR-Tools proporciona una herramienta principal para resolver este tipo de problemas de optimización:

El problema

En el artículo de introducción (programación lineal), abordamos un caso descrito en el libro Applied Mathematical Programming, de Bradley, Hax, and Magnanti (Addison-Wesley, 1977), del MIT (Cápitulo 2 página 50). Con fines prácticos, hemos adaptado dicho problema, incorporando nuevas restricciones y modificando algunos datos del modelo original, con el propósito de evidenciar la diferencia en los resultados obtenidos por medio del uso de variables continuas, enteras y mixtas.

El propietario de una tienda que produce remolques para automóviles desea determinar la mejor combinación para sus tres productos: remolques de plataforma plana, remolques económicos y remolques de lujo. Su taller se limita a trabajar 24 días al mes en el trabajo de los metales y 60 días al mes en el trabajo de la madera para estos productos.

Existe un contrato vigente, mediante el cual, el propietario deberá entregar como mínimo 4 remolques tipo económico cada mes.

La siguiente tabla indica los datos de producción de los remolques.

Modelamiento del problema

Sean las variables de decisión del problema:

x₀ = Número de remolques de plataforma plana producidos por mes

x₁ = Número de remolques económicos producidos por mes

x₂ = Número de remolques de lujo producidos por mes

Suponiendo que los costos de la capacidad de trabajo en metal y madera sean fijos, el problema se convierte en un problema de maximización:

Zmax = 6x₀ + 14x₁ + 13x₂

Sujeto a las siguientes restricciones de capacidad:

0,6x₀ + 2x₁ + x₂ <= 24,

x₀ + 2x₁ + 4x₂ <= 60,

Sujeto a la siguiente restricción de demanda mínima (contrato)

x₁ >= 4

Sujeto a las siguientes restricciones de no-negatividad:

x₀ >= 0,

x₁ >= 0,

x₂ >= 0,

Podemos, del mismo modo, establecer un par de variables que correspondan a las horas ociosas para las dos tareas establecidas (metal y madera):

x₃ = Número de horas ociosas en el trabajo en metal al mes,

x₄ = Número de horas ociosas en el trabajo en madera al mes,

Reescribimos las restricciones (adicionando las variables de horas ociosas). Podemos observar que las inecuaciones ahora serán igualdades, para que de esta forma ahora podamos tener información relacionada a los recursos. En otras palabras, lo que se utiliza (horas productivas) + lo que sobre (horas ociosas) = tiempo disponible:

0,6x₀ + 2x₁ + x₂ + x₃ = 24,

x₀ + 2x₁ + 4x₂ + x₄ = 60,

x₀ >= 0,

x₁ >= 0,

x₂ >= 0,

x₃ >= 0,

x₄ >= 0,

Así entonces, tenemos el problema completamente modelado.

Resolución del modelo mediante Google OR-Tools

De acuerdo a lo mencionado en el artículo de introducción a Google OR-Tools, esta herramienta soporta múltiples lenguajes de programación, en esta ocasión, haremos uso del lenguaje de programación Python.

Paso 1: Importar la librería

El siguiente fragmento de código importa las librerías necesarias:

# Importar la librería de Google OR-Tools
from ortools.linear_solver import pywraplp

Paso 2: Declarar el solucionador

El siguiente fragmento de código declara el solucionador SCIP (Solving Constraint Integer Programs), un solucionador de código abierto disponible que permite resolver problemas lineales mixtos (Google OR-Tools posee múltiples solucionadores):

solver = pywraplp.Solver.CreateSolver('SCIP')

Paso 3: Crear las variables del modelo

El siguiente fragmento de código permite crear las variables del modelo. La sintaxis permite declarar la naturaleza de cada una de las variables y el rango de valores permitidos (restricciones de no-negatividad).

    x0 = solver.IntVar(0, solver.infinity(), 'x0')
    x1 = solver.IntVar(0, solver.infinity(), 'x1')
    x2 = solver.IntVar(0, solver.infinity(), 'x2')
    x3 = solver.NumVar(0, solver.infinity(), 'x3')
    x4 = solver.NumVar(0, solver.infinity(), 'x4')

A partir de la declaración de estas variables, el modelo corresponde a un problema de programación lineal mixta.

Paso 4: Definir las restricciones del modelo

El siguiente fragmento de código define las restricciones del modelo:

    # Restricción 0: 0.6x0 + 2x1 + x2 + x3 = 24. (horas metales)
    solver.Add(0.6 * x0 + 2 * x1 + x2 + x3 == 24.0)

    # Restricción 1: x0 + 2x1 + 4x2 + x4 = 60. (horas madera)
    solver.Add(x0 + 2 * x1 + 4 *x2 + x4 == 60.0)

    # Restricción 3: x1 >= 4. (demanda mínima)
    solver.Add(x1 >= 4.0)

Paso 5: Definir la función objetivo del modelo

El siguiente fragmento de código define la función objetivo del modelo:

    # Función objetivo (max): 6x0 + 14x1 + 13x2
    solver.Maximize(6 * x0 + 14 * x1 + 13 * x2)

Paso 6: Invocar el solucionador

El siguiente fragmento de código sirve para invocar el solucionador del modelo:

status = solver.Solve()

Paso 7: Definir las salidas del solucionador

    if status == pywraplp.Solver.OPTIMAL:
        print('Solución:')
        print('Valor objetivo =', solver.Objective().Value())
        print('Número de remolques de plataforma plana producidos por mes =', x0.solution_value())
        print('Número de remolques económicos producidos por mes =', x1.solution_value())
        print('Número de remolques de lujo producidos por mes =', x2.solution_value())
        print('Número de horas ociosas en el trabajo en metal al mes =', x3.solution_value())
        print('Número de horas ociosas en el trabajo en madera al mes =', x4.solution_value())
    else:
      if status == solver.FEASIBLE:
        print('Se encontró una solución potencialmente subóptima.')
      else:
        print('El problema no tiene solución óptima.')

    # Información avanzada del solucionador

    print('\nUso avanzado:')
    print('Problema resuelto en %f milisegundos' % solver.wall_time())
    print('Problema resuelto en %d iteraciones' % solver.iterations())

Es posible que el desarrollo de los siete pasos anteriores demande algún grado de complejidad subyacente del uso de un lenguaje de programación; sin embargo, es preciso mencionar que, el modelo anterior queda perfectamente configurado, y puede replicarse con modificaciones menores en múltiples problemas de optimización lineal.

¿Cómo ejecutar el modelo?

Alternativa 1, ejecución en nuestro equipo:

Lo primero que debemos considerar, en el caso de que queramos ejecutar este código en nuestro equipo, es que es preciso contar con la instalación de Python en nuestro equipo de cómputo, así mismo debemos contar con la última versión del comando pip y por supuesto, el software OR-Tools. Una guía detallada de la instalación de estos requerimientos la podrás encontrar en el siguiente enlace:

Instalación de OR-Tools para Python

Ahora, lo recomendable es trabajar con algún editor de código práctico (IDE), por ejemplo: Sublime Text, o Spyder (Una herramienta más completa y por ende más robusta y pesada).

Alternativa 2, ejecución en un entorno virtual (Recomendado):

Podemos utilizar del mismo modo, un entorno virtual. En este caso recomendamos el uso de Colaboratory de Google, un entorno que cuenta con todas las herramientas necesarias para nuestros desarrollos. No tendremos que instalar nada en nuestro equipo, y aprovecharemos la potencia de las máquinas de Google.

El código completo de nuestro desarrollo lo presentamos a continuación. También puedes ver el cuaderno de este módulo en nuestro Colaboratory: Programación Lineal Mixta.

#Adaptado de: Bradley, Hax, and Magnanti, 'Applied Mathematical Programming', Chapter 2
#Nuevo caso y modelo: Salazar, ingenieriaindustrialonline.com (Programación lineal mixta)

# Importar la librería de Google OR-Tools
from ortools.linear_solver import pywraplp

def LinearProgrammingExample():
    solver = pywraplp.Solver.CreateSolver('SCIP')

    x0 = solver.IntVar(0, solver.infinity(), 'x0')
    x1 = solver.IntVar(0, solver.infinity(), 'x1')
    x2 = solver.IntVar(0, solver.infinity(), 'x2')
    x3 = solver.NumVar(0, solver.infinity(), 'x3')
    x4 = solver.NumVar(0, solver.infinity(), 'x4')

    print('Número de variables =', solver.NumVariables())

    # Restricción 0: 0.5x0 + 2x1 + x2 + x3 = 24.
    solver.Add(0.6 * x0 + 2 * x1 + x2 + x3 == 24.0)

    # Restricción 1: x0 + 2x1 + 4x2 + x4 = 60.
    solver.Add(x0 + 2 * x1 + 4 *x2 + x4 == 60.0)

    # Restricción 3: x1 >= 4.
    solver.Add(x1 >= 4.0)

    print('Número de restricciones =', solver.NumConstraints())

    # Función objetivo (max): 6x0 + 14x1 + 13x2
    solver.Maximize(6 * x0 + 14 * x1 + 13 * x2)

    # Declarar el solucionador.
    status = solver.Solve()

    # Declarar las salidas del solucionador
    if status == pywraplp.Solver.OPTIMAL:
        print('Solución:')
        print('Valor objetivo =', solver.Objective().Value())
        print('Número de remolques de plataforma plana producidos por mes =', x0.solution_value())
        print('Número de remolques económicos producidos por mes =', x1.solution_value())
        print('Número de remolques de lujo producidos por mes =', x2.solution_value())
        print('Número de horas ociosas en el trabajo en metal al mes =', x3.solution_value())
        print('Número de horas ociosas en el trabajo en madera al mes =', x4.solution_value())
    else:
      if status == solver.FEASIBLE:
        print('Se encontró una solución potencialmente subóptima.')
      else:
        print('El problema no tiene solución óptima.')

    # Información avanzada del solucionador

    print('\nUso avanzado:')
    print('Problema resuelto en %f milisegundos' % solver.wall_time())
    print('Problema resuelto en %d iteraciones' % solver.iterations())

LinearProgrammingExample()

Al ejecutar nuestro desarrollo en Colaboratory, tenemos:

Hemos modificado la naturaleza de las variables con el objetivo de mostrar los resultados a partir de tres escenarios: variables mixtas, enteras y continuas:

El anterior, es un problema lineal que representa un caso sencillo, sin embargo, los resultados obtenidos presentan variaciones menores de acuerdo a la naturaleza de las variables consideradas.

Sin embargo, en modelos robustos, estas variaciones pueden ser considerables y determinantes. Por tal razón, la programación lineal mixta ofrece la posibilidad de modelar variables cuya naturaleza refleje con precisión la realidad que se pretende representar.

Ahora bien, el modelo de optimización lineal mixta y el script del solucionador quedaron desarrollados en un lenguaje de programación estándar y ampliamente utilizado. Desde luego, las posibilidades de integrar datos de entrada y procesar los datos de salidas son interesantes. Por ejemplo, es posible desarrollar un script mediante el cual el código ya desarrollado tome los datos de entrada desde un archivo de Excel, o desde un servidor externo.

También, es posible desarrollar una interfaz amigable desde la cual se ingrese la información; o vincular los datos de salida con algún modelo o documento determinado.