Otimizando com Pyomo: Guia passo a passo completo

Configurando o Pyomo

Antes de mergulharmos na modelagem, precisamos configurar nosso ambiente instalando o Pyomo e escolhendo um solucionador apropriado.

1. Pré-requisitos

Para usar o Pyomo, você deve ter Python 3.6 ou superior. O Pyomo pode ser instalado através do pip.

pip install pyomo

Este tutorial foi criado usando a versão 6.8.0 do Pyomo.

import pyomo
print(pyomo.__version__)

Saída:

>>> 6.8.0

2. Escolhendo e instalando o solucionador correto

Na otimização, os solucionadores são essenciais pois eles são os algoritmos que encontram a solução ótima para o problema que você definiu. Diferentes solucionadores são melhores adaptados dependendo do tipo de problema (por exemplo, linear, não linear, inteiro). Pyomo é uma ferramenta de modelagem que depende de solucionadores externos para realizar as computações reais.

Vamos revisitar alguns dos solucionadores mais comuns.

Solucionadores de código aberto

1. GLPK (GNU Linear Programming Kit)

GLPK é uma ferramenta popular para resolver problemas de programação linear (LP) e de programação linear mista com inteiros (MIP).

É uma excelente escolha para tarefas básicas de otimização linear e é amplamente usada em aplicações acadêmicas e industriais.

Instalação

• Windows: Faça download e instale GLPK.
• macOS: 

brew install glpk

• Linux: 

sudo apt-get install glpk-utils

2. CBC (Coin-or Branch and Cut)

CBC é um solucionador de código aberto para problemas de programação linear (LP) e de programação inteira mista (MIP).

Ele oferece funcionalidades avançadas e melhores desempenhos, em alguns casos, em comparação com o GLPK, tornando-se uma boa opção para tarefas de otimização mais complexas.

O CBC pode ser instalado através do gerenciador de pacotes conda.

conda install -c conda-forge coincbc

3. IPOPT (Interior Point OPTimizer)

IPOPT é um poderoso solucionador projetado para problemas de programação não linear (NLP) em grande escala. 

É particularmente adequado para lidar com modelos complexos e não lineares, tornando-o uma excelente escolha para problemas além da optimização linear.

IPOPT pode ser instalado via o gerenciador de pacotes conda.

!conda install -c conda-forge ipopt

Solucionadores comerciais

1. CPLEX

CPLEX é um solucionador de otimização de ponta de gênero que gerencia eficientemente problemas de programação linear (LP), programação inteira mista (MIP) e problemas de programação quadrática (QP).

É necessário um licença da IBM, mas está disponível de graça para usuários acadêmicos, tornando-se uma excelente escolha para fins de pesquisa e educacionais.

2. Gurobi

Gurobi é um solucionador comercial líder conhecido por sua velocidade e eficiência em resolver problemas LP, MIP, QP e programação não linear (NLP).

Como o CPLEX, ele exige uma licença mas oferece acesso grátis para usuários acadêmicos. Portanto, é uma ferramenta padrão da indústria para otimização avançada.

Solucionadores open-source vs. comerciais

Solucionadores de código aberto, como o GLPK e o CBC, são gratuitos e suficientes para a maioria das necessidades básicas de otimização. Eles são excelentes escolhas para projetos de pequena escala e fins educacionais.

No entanto, solucionadores comerciais, como o CPLEX e o Gurobi, oferecem normalmente melhor desempenho, especialmente para problemas maiores e mais complexos. Estes solucionadores incluem funcionalidades avançadas, incluindo suporte aprimorado para programação quadrática e não linear, e são otimizados para aplicações industriais de grande escala.

Enquanto solucionadores de código aberto podem lidar com muitas tarefas de otimização rotineira, solucionadores comerciais são frequentemente uma melhor escolha quando se trata de requisitos mais complexos e de alto desempenho.

Note que os solucionadores comerciais exigem uma licença, embora sejam disponíveis de graça para usuários acadêmicos.

Agora, vejamos como configurar um solucionador em Pyomo. Na sequência, eu usarei o GLPK neste caso.

3. Configurando um solucionador em Pyomo

Primeiro, certifique-se que o executável do solucionador esteja em seu caminho de execução do sistema após a instalação.

Em seguida, crie um script Python e adicione o seguinte:

from pyomo.environ import SolverFactory
solver = SolverFactory('glpk')

Para confirmar que ambos Pyomo e o seu solucionador estão instalados corretamente, vamos resolver um problema de teste simples.

Problema de teste: programa linear simples.

Objetivo: Minimizar Z=x+y

Subjeto a:

• x + 2y ≥ 4
• x – y ≤ 1
• x ≥ 0
• y ≥ 0

Este problema é sobre encontrar o menor valor possível de Z, que é a soma de duas variáveis, x e y. No entanto, x e y devem atender a certas condições.

Primeiro, quando você adiciona x e duas vezes y, o resultado deve ser pelo menos 4. Segundo, x menos y deve ser menor ou igual a 1. Finalmente, ambos x e y devem ser zero ou números positivos (eles não podem ser negativos). 

O objetivo é encontrar valores de x e y que satisfazem estas condições, enquanto torna Z o menor possível.

Implementando usando Pyomo:

import pyomo.environ as pyo
# Criar um modelo
model = pyo.ConcreteModel()
# Definir variáveis
model.x = pyo.Var(within=pyo.NonNegativeReals)
model.y = pyo.Var(within=pyo.NonNegativeReals)
# Definir objetivo
model.obj = pyo.Objective(expr=model.x + model.y, sense=pyo.minimize)
# Definir restrições
model.con1 = pyo.Constraint(expr=model.x + 2 * model.y >= 4)
model.con2 = pyo.Constraint(expr=model.x - model.y <= 1)
# Selecionar solucionador
solver = pyo.SolverFactory('glpk')
# Solucionar o problema
result = solver.solve(model)
# Exibir resultados
print('Status:', result.solver.status)
print('Termination Condition:', result.solver.termination_condition)
print('Optimal x:', pyo.value(model.x))
print('Optimal y:', pyo.value(model.y))
print('Optimal Objective:', pyo.value(model.obj))

Se tudo estiver funcionando corretamente, a saída esperada será:

Status: ok
Termination Condition: optimal
Optimal x: 0.0
Optimal y: 2.0
Optimal Objective: 2.0

Vamos passar pelaqueles códigos acima: Primeiro, ele define duas variáveis, x e y, que podem assumir apenas valores não negativos. O objetivo do modelo é minimizar a soma de x e y (x + y). O código define o solucionador como glpk para encontrar os valores ótimos de x e y que satisfazem essas restrições e minimizam o objetivo.

Após executar o código, descobrimos que os valores ótimos para as variáveis são x = 0.0 e y = 2.0, que minimizam a função objetivo Z = x + y. Portanto, o valor mínimo da função objetivo é 2.0, que satisfaz as restrições dadas.

Modelagem básica com Pyomo

Entender como definir os componentes básicos de um modelo de otimização em Pyomo é necessário para configurar e resolver problemas de otimização com eficiência.

1. Definição de variáveis

Variáveis representam as decisões que precisam ser feitas em um problema de otimização. Em Pyomo, variáveis são as quantidades que o solucionador irá ajustar para otimizar a função objetivo enquanto satisfaz todas as restrições.

Variáveis escalares

Uma variável escalar é uma variável simples que não é indexada sobre qualquer conjunto. Para definir uma variável escalar em Pyomo, você usa a classe Var do módulo pyomo.environ.

from pyomo.environ import Var
model.x = Var()

Primeiro, importamos a Var e criamos uma variável x usando Var(). Esta variável não tem limites especificados, o que significa que ela pode assumir qualquer valor real, a menos que seja restringida de outra forma no modelo.

Adicionando limites

Você pode restringir os valores que uma variável pode assumir especificando limites. Os limites são definidos como uma tupla (lower_bound, upper_bound):

from pyomo.environ import Var
model.x = Var(bounds=(0, None))

Especificando domínios

O Pyomo fornece domínios predefinidos que você pode usar para especificar o tipo de valores que uma variável pode assumir, como NonNegativeReals, Integers ou Binary:

from pyomo.environ import Var, NonNegativeReals
model.x = Var(domain=NonNegativeReals)

Variáveis indexadas

Quando se trata de múltiplas variáveis que são semelhantes por natureza, como variáveis que representam períodos de tempo diferentes ou itens, é eficiente usar variáveis indexadas. Variáveis indexadas são variáveis definidas sobre um conjunto.

import pyomo.environ as pyo
model.I = pyo.Set(initialize=[1, 2, 3])
model.y = pyo.Var(model.I, domain=pyo.NonNegativeReals)

Suponha que você esteja modelando as quantidades de produção de três produtos. Você pode definir:

model.Products = pyo.Set(initialize=['A', 'B', 'C'])
model.production = pyo.Var(model.Products, domain=pyo.NonNegativeReals)

Agora, model.production['A'], model.production['B'], e model.production['C'] representam as quantidades de produção para os produtos A, B e C, respectivamente.

2. Definição de objetivos

A função objetivo é o que estamos tentando otimizar (minimizar ou maximizar). Ela define o objetivo do modelo, como minimizar custos ou maximizar lucros, e é normalmente expressa como uma equação matemática envolvendo as variáveis de decisão.

Elas são definidas usando a classe Objective:

from pyomo.environ import ConcreteModel, Var, Objective, minimize, maximize, NonNegativeReals
# Criar um modelo
model = ConcreteModel()
# Definir variáveis
model.x = Var(within=NonNegativeReals)
model.y = Var(within=NonNegativeReals)
# Minimização (custo)
model.cost = Objective(expr=2 * model.x + 3 * model.y, sense=minimize)
# Quando a maximização de lucros - (pode haver apenas um objetivo de cada vez)
# model.profit = Objective(expr=5 * model.x + 4 * model.y, sense=maximize)

3. Adicionando restrições

As restrições definem as limitações ou requisitos do problema:

from pyomo.environ import Constraint
model.con1 = Constraint(expr=model.x + model.y >= 10)

O exemplo acima define uma restrição no modelo Pyomo usando a classe Constraint. A restrição model.con1 especifica que a soma das variáveis x e y deve ser maior ou igual a 10.

4. Parametrizando modelos

Parâmetros são valores fixos usados no modelo para representar quantidades conhecidas ou constantes que não mudam durante o processo de otimização.

Eles ajudam a definir as relações entre variáveis e restrições, fornecendo estrutura ao modelo através da incorporação de dados reais ou pressupostos:

from pyomo.environ import Param
model.p = Param(initialize=5)

O código acima define um parâmetro p no modelo Pyomo usando a classe Param e inicializando-o com um valor fixo de5. O parâmetro p pode agora ser usado no modelo para representar um valor constante que não muda durante o processo de otimização.

Agora, vamos trabalhar em um problema de otimização de ponta a ponta!

Exemplo de Custo de Ponta a Ponta com Pyomo

Vamos ver um exemplo de resolução de problema de otimização de ponta a ponta usando Pyomo. Vamos modelar uma situação real onde uma fábrica produz dois produtos e o objetivo é maximizar o lucro considerando as restrições de tempo de máquina.

1. Declaração do problema

Uma fábrica produz dois produtos, P1 e P2. O lucro por unidade é:

• P1: $40
• P2: $50

Tempo de máquina disponível:

• Máquina A: 100 horas
• Máquina B: 80 horas
• Máquina C: 90 horas

Tempo necessário por unidade:

Objetivo: Maximizar lucro.

Variáveis de decisão:

• x₁: Unidades de P1 a produzir.
• x₂: Unidades de P2 para produção.

2. Formulação matemática

Função objetivo:

Maximizar Z = 40x₁ + 50x₂

Restrições:

1. Capacidade da máquina A: 1x₁ + 2x₂ ≤ 100
2. Capacidade da máquina B: 2x₁ + 1x₂ ≤ 80
3. Capacidade da máquina C: 3x₂ ≤ 90
4. Non-negativity: x₁, x₂ ≥ 0

3. Implementação

Com base no objetivo e nas restrições do problema, aqui está o código Python para modelá-lo, novamente, usando GLPK.

# Passo 1: Importar Bibliotecas
import pyomo.environ as pyo
# Passo 2: Criar um Modelo Concreto
model = pyo.ConcreteModel()
# Passo 3: Definir Variáveis de Decisão (Unidades de P1 e P2 para produzir)
model.x1 = pyo.Var(within=pyo.NonNegativeReals)
model.x2 = pyo.Var(within=pyo.NonNegativeReals)
# Passo 4: Definir a Função Objetivo (Maximizar lucro)
model.profit = pyo.Objective(expr=40 * model.x1 + 50 * model.x2, sense=pyo.maximize)
# Passo 5: Definir Constraintes
# Constraintes de capacidade da máquina A: 1x1 + 2x2 <= 100
model.machine_a = pyo.Constraint(expr=1 * model.x1 + 2 * model.x2 <= 100)
# Constraintes de capacidade da máquina B: 2x1 + 1x2 <= 80
model.machine_b = pyo.Constraint(expr=2 * model.x1 + 1 * model.x2 <= 80)
# Constraintes de capacidade da máquina C: 3x2 <= 90
model.machine_c = pyo.Constraint(expr=3 * model.x2 <= 90)
# Passo 6: Solucionar o Modelo usando o solucionador GLPK
solver = pyo.SolverFactory('glpk')
result = solver.solve(model)
# Passo 7: Analisar Resultados
# Exibir Status do Solucionador e Condição de Encerramento
print('Solver Status:', result.solver.status)
print('Termination Condition:', result.solver.termination_condition)
# Obter e exibir os valores ótimos para x1, x2, e o lucro máximo
x1_opt = pyo.value(model.x1)
x2_opt = pyo.value(model.x2)
profit_opt = pyo.value(model.profit)
print(f'Optimal production of P1 (x1): {x1_opt}')
print(f'Optimal production of P2 (x2): {x2_opt}')
print(f'Maximum Profit: ${profit_opt}')

Output:

>>> Solver Status: ok
>>> Termination Condition: optimal
>>> Optimal production of P1 (x1): 25.0
>>> Optimal production of P2 (x2): 30.0
>>> Maximum Profit: $2500.0

No código acima, definimos um modelo de otimização linear para maximizar o lucro da produção de dois produtos (P1 e P2). A função objetivo é definida para maximizar o lucro, com cada unidade de P1 a contribuir com $40 e cada unidade de P2 a contribuir com $50.

Impõemos três restrições que representam os limites de tempo de máquina para as Máquinas A, B, e C.

Finalmente, usamos o solucionador GLPK para resolver o problema.

A resposta final é a de produzir 25 unidades de P1 e 30 unidades de P2 onde o nosso lucro máximo será de $2,500.

Recursos Avançados em Pyomo

Na seção anterior, vimos como é fácil implementar um problema de otimização de ponta a ponta com Pyomo. No entanto, a maioria dos problemas da vida real não é fácil de resolver.

Nesta seção, apresento algumas funcionalidades avançadas que você pode usar para resolver cenários mais complexos.

1. Otimização Não Linear

Otimização Não Linear minimiza ou maximiza uma função objetivo não linear, sujeita a restrições não lineares. Vamos olhar para um exemplo onde minimizamos a soma das diferenças quadradas sujeitas a uma restrição circular.

Enunciado do Problema

Minimizar o objetivo:Z = (x – 1)² + (y – 2)²

Subject to:

• x² + y² ≤ 4
• x, y ≥ 0

Em Pyomo, podemos definir as variáveis de decisão x e y com limites de 0 para garantir não negatividade. A função objetivo é escrita como a soma das diferenças quadradas de pontos específicos, e a restrição garante que a solução esteja dentro de um círculo de raio 2.

Neste caso, o solucionador IPOPT é adequado para sua capacidade de resolver problemas de otimização não linear:

import pyomo.environ as pyo
model = pyo.ConcreteModel()
# Defina variáveis com limites inferiores
model.x = pyo.Var(bounds=(0, None))
model.y = pyo.Var(bounds=(0, None))
# Função objetivo: minimizar (x - 1)² + (y - 2)²
model.obj = pyo.Objective(expr=(model.x - 1)**2 + (model.y - 2)**2, sense=pyo.minimize)
# Restrição: x² + y² ≤ 4 (circulo de raio 2)
model.circle = pyo.Constraint(expr=model.x**2 + model.y**2 <= 4)
solver = pyo.SolverFactory('ipopt')
result = solver.solve(model)
print('Optimal x:', pyo.value(model.x))
print('Optimal y:', pyo.value(model.y))
print('Minimum Z:', pyo.value(model.obj))

2. Programação mista com inteiros (MIP)

A programação mista com inteiros é usada quando algumas variáveis de decisão são inteiros ( frequentemente binários) enquanto outras são contínuas. É valioso para problemas de decisão como localização de instalações e planejamento de produção.

Declaração do problema

Uma empresa deve decidir se deve abrir armazéns em locais A, B, e C. O objetivo é minimizar o custo total, que inclui os custos fixes de abertura de armazéns e os custos de transporte.

Começamos inicializando os dados, incluindo os custos fixes de abertura de armazéns, custos de transporte, limites de capacidade e a demanda total:

locations = ['A', 'B', 'C']
FixedCost = {'A': 1000, 'B': 1200, 'C': 1500}
TransportCost = {'A': 5, 'B': 4, 'C': 6}
Capacity = {'A': 100, 'B': 80, 'C': 90}
Demand = 150
model = pyo.ConcreteModel()
# Variável binária: 1 se o armazém estiver aberto, 0 caso contrário
model.y = pyo.Var(locations, domain=pyo.Binary)
# Variável contínua: quantidade de mercadorias transportadas
model.x = pyo.Var(locations, domain=pyo.NonNegativeReals)
model.cost = pyo.Objective(
    expr=sum(FixedCost[i] * model.y[i] + TransportCost[i] * model.x[i] for i in locations),
    sense=pyo.minimize
)
# Restrição de demanda
model.demand = pyo.Constraint(expr=sum(model.x[i] for i in locations) >= Demand)
# Restrições de capacidade
def capacity_rule(model, i):
    return model.x[i] <= Capacity[i] * model.y[i]
model.capacity = pyo.Constraint(locations, rule=capacity_rule)
solver = pyo.SolverFactory('cbc')
result = solver.solve(model)
for i in locations:
    print(f"Warehouse {i}: Open={pyo.value(model.y[i])}, Transported={pyo.value(model.x[i])}")
print('Minimum Total Cost:', pyo.value(model.cost))

O modelo inclui dois tipos de variáveis de decisão: uma variável binária y que representa se um armazém está aberto (1 se aberto, 0 caso contrário) e uma variável contínua x que representa a quantidade de mercadorias transportadas de cada armazém.

A função objetivo soma os custos fixos e de transporte de cada armazém e minimiza o total. As restrições garantem que o total de mercadorias transportadas atende à demanda e que a capacidade de cada armazém não é excedida se ele for aberto.

3. Manuseio de múltiplos objetivos

Às vezes, problemas de otimização envolvem múltiplos objetivos que podem entrar em conflito, como maximizar o lucro enquanto minimiza o impacto ambiental. Uma abordagem comum é o método de soma ponderada, onde cada objetivo é atribuído um peso para equilibrar sua importância.

Enunciado do problema

Nosso objetivo é maximizar o lucro enquanto minimiza o impacto ambiental:

• Lucro:Z₁ = 3x + 5y
• Impacto ambiental:Z₂ = 2x + y

Podemos combinar estes objetivos usando as pesos w1=0.6, w2=0.4, onde o objetivo total se torna uma soma ponderada:

w1 = 0.6
w2 = 0.4
model.obj = pyo.Objective(
    expr=w1 * (3 * model.x + 5 * model.y) - w2 * (2 * model.x + model.y),
    sense=pyo.maximize
)

Neste objetivo combinado, nós maximizamos o lucro enquanto minimizamos o impacto ambiental alterando as pesagens.

4. Utilizar fontes de dados externas

Ao trabalhar com grandes conjuntos de dados, importar dados de fontes externas, como arquivos CSV, é frequentemente útil. O Pyomo funciona bem com o Pandas para ler e usar dados externos.

Nós podemos ler um arquivo CSV usando o Pandas e usar os dados para inicializar conjuntos e parâmetros no nosso modelo:

import pandas as pd
data = pd.read_csv('parameters.csv')
# Definir conjunto a partir de dados CSV
model.I = pyo.Set(initialize=data['index'].unique())
# Definir parâmetro inicializado a partir de dados CSV
param_dict = data.set_index('index')['value'].to_dict()
model.param = pyo.Param(model.I, initialize=param_dict)

Dicas e Melhores Práticas para Usar Pyomo

Ao trabalhar com Pyomo, manter seus modelos eficientes, bem documentados e fáceis de diagnosticar é importante.

1. Depuração e resolução de problemas

Enquanto constrói modelos de otimização em Pyomo, é comum encontrar problemas como soluções infeasíveis, falhas do solucionador ou resultados incorretos. Aqui estão algumas melhores práticas para depuração:

• Verificação de restrições: Confera suas restrições se o modelo não estiver produzindo uma solução viável. Restrições rígidas podem tornar um problema inviável. Use o método .display() do Pyomo para imprimir os valores das variáveis e restrições para verificar que eles estão se comportando como esperado.
• Saída do solucionador: Ative logs detalhados do solucionador passando tee=True quando chamar o método solve(). Isto pode fornecer insights sobre onde o solucionador poderá ter dificuldades, como variáveis desvinculadas ou infeasibilidade.
• Teste primeiro modelos simples: Quando se trata de modelos complexos, teste uma versão simplificada. Isso pode ajudar a isolar potenciais problemas sem a sobrecarga de um modelo totalmente especificado.

A depuração de problemas é muito mais fácil se você abordá-la de forma sistemática, analisando restrições, a função objetivo e a resposta do solucionador.

2. Eficiência no modelagem

Problemas de otimização podem se tornar computacionalmente caros conforme o tamanho do modelo aumenta. Para garantir uma modelagem eficiente, considere as seguintes dicas:

• Usar esparsidade: Evite percorrer índices desnecessários ao definir restrições ou objetivos. Aproveitar a esparsidade no seu problema reduz o tempo de cálculo.
• Variáveis binárias vs. contínuas: Onde possível, reduza o número de variáveis binárias ou inteiras. Variáveis contínuas são mais fáceis para os solucionadores lidarem, resultando em soluções mais rápidas.
• Formulação de restrições: Mantenha as restrições o mais simples possível, tanto na forma matemática quanto na implementação. Evite não necessidades não lineares e descomponha restrições complexas em pequenas e gerenciáveis.

Modelos eficientes soluções mais rápido e são mais fáceis de depurar e manter.

3. Documentação e manutenção

Manter modelos Pyomo bem documentados é uma boa prática para uso longo e colaboração. Boa documentação torna também mais fácil voltar e atualizar modelos ao longo do tempo:

• Usar comentários em linha: Sempre adicionar comentários para explicar o propósito de variáveis, restrições e função objetivo. Isto é especialmente importante em modelos de otimização onde a lógica pode não ser imediatamente óbvia.
• Modularize seu código: Separe seu modelo em seções lógicas ou até mesmo em funções separadas. Esse approach modular pode melhorar a legibilidade e tornar mais fácil a depuração e a modificação de partes específicas do modelo.
• Rastreie as mudanças no modelo: Mantenha um histórico de versões do seu modelo, especialmente se ele estiver evoluindo. Use ferramentas de controle de versão, como Git, para rastrear as mudanças e garantir que quaisquer atualizações ou melhorias possam ser voltadas.

Documentação apropriada e código estruturado tornarão seus modelos Pyomo mais acessíveis a colaboradores futuros e tornará mais fácil escalar ou modificar conforme suas exigências evoluem.

Conclusão

Pyomo é uma ferramenta poderosa e flexível para construir e solucionar modelos de otimização em Python. Durante este tutorial, explorámos como o Pyomo permite que usuários modelem vários problemas de otimização, desde programação linear até programação não linear e mista com inteiros.

Com sua sintaxe amigável e integração com solucionadores, o Pyomo torna a formulação e solução de problemas reais de otimização acessíveis tanto para iniciantes quanto para usuários avançados.

Se você está interessado em aprender mais sobre a solução de problemas reais com otimização, verifique o curso gratuito Introdução à Otimização em Python no DataCamp!