Programação e sequenciamento das frentes de colheita de cana-de-açúcar: modelo e métodos de solução para problemas de grande porte

Resumo:

Em um estudo anterior recente, modelos de otimização para a programação e o sequenciamento das frentes de colheita de cana-de-açúcar foram propostos. Esse é um problema agrícola e logístico complexo que envolve vários fatores, tais como o estágio de maturação da matéria-prima, a colheita na unidade agrícola e o transporte dessa matéria-prima para a unidade industrial, bem como a capacidade de moagem da unidade industrial. No presente estudo, aplica-se um dos modelos de otimização do estudo anterior para representar esse problema por meio de um modelo de programação inteira mista (PIM) de dimensionamento e sequenciamento de lotes da produção em máquinas paralelas com custos e tempos de setup dependentes da sequência. Propõem-se métodos baseados em heurísticas PIM para resolver esse modelo em uma situação real de uma safra de uma empresa típica do setor, inspirados em uma heurística de agregação de blocos de colheita, em heurísticas de construção do tipo relax-and-fix e heurísticas de melhoria do tipo fix-and-optimize. Para comparar os desempenhos desses métodos heurísticos foram realizados vários experimentos com diferentes combinações e variações desses métodos, e três abordagens foram capazes de gerar soluções de boa qualidade, sendo que uma delas é aqui detalhada e analisada, com resultados promissores para apoiar decisões de programação e sequenciamento das frentes de colheita de cana-de-açúcar.

Palavras-chave:
Programação de colheita de cana-de-açúcar; Dimensionamento e sequenciamento de lotes de produção; Programação inteira mista; Heurísticas relax-and-fix e fix-and-optimize

Abstract:

In a recent study, optimization models were proposed for programming and scheduling sugarcane harvesting fronts. This is a complex agricultural and logistic problem comprising various factors, such as raw material maturation stage, harvesting at the agricultural unit, transporting of raw material to the plant, and milling capacities of the plant. In this study, one of the optimization models previously studied was used to represent this problem using Mixed Integer Programming (MIP) of a lot sizing and scheduling model in parallel machines with sequence dependent setup times and costs. The proposed methods are based on MIP heuristics to solve this model in a real situation of a harvest season of a typical company from this sector inspired by harvest block aggregation heuristics, relax-and-fix constructive heuristics, and fix-and-optimize improvement heuristics. To compare the performance of the heuristic methods, various experiments were conducted using different combinations and variations of these methods. Three approaches were able to produce good quality solutions. One of them is described in detail and analyzed in this study, showing promising results in terms of making programming and scheduling decisions concerning sugarcane harvesting fronts.

Keywords:
Sugarcane harvest programming; Production lot sizing and scheduling; Mixed integer programming; Relax-and-fix and fix-and-optimize heuristics

1 Introdução

Abordagens de pesquisa operacional aplicadas no planejamento da produção e da logística de agroindústrias vêm sendo largamente estudadas na literatura (Ahumada & Villalobos, 2009; Junqueira & Morabito, 2012; Plà et al., 2014; Plà-Aragones, 2015). Nesse contexto, o planejamento de colheita de safra de cana-de-açúcar é um problema complexo que envolve vários fatores, tais como o estágio de maturação da matéria-prima, sua colheita na unidade agrícola e seu transporte para a unidade industrial, bem como a capacidade de moagem da unidade industrial, sendo ele fundamental para a integração das etapas agrícola e industrial. Conforme Junqueira & Morabito (2017), parte dessa complexidade é fruto da divisão da operação agrícola em diversas frentes de colheita e das dificuldades de programação e sequenciamento dessas frentes ao longo do tempo, de maneira a manter uma capacidade média de transporte de cana para a usina, conforme as frentes vão mudando de blocos de colheita. A manutenção dessa capacidade média de transporte poderia se tornar inviável caso todo o esforço de colheita estivesse concentrado em apenas poucos pontos (por exemplo, todas as frentes localizadas num único bloco de colheita), pois ora os recursos poderiam estar em excesso, se a maior parte dessas frentes estivesse localizada muito próxima da usina, ora poderiam estar em falta, se essas frentes estivessem muito distantes.

Vários trabalhos trataram do planejamento de colheita de safra de cana-de-açúcar, tais como Higgins et al. (1998), Higgins (1999, 2002), Higgins & Muchow (2003) e Higgins et al. (2004a, b), no contexto australiano; Grunow et al. (2007), no venezuelano, e Jena & Poggi (2013), no brasileiro. Entretanto, nesses estudos, a divisão em frentes de colheita ou foi desconsiderada, ou predefinida por bloco de colheita. Em um estudo anterior, Junqueira & Morabito (2017) se inspiraram no Problema Geral de Dimensionamento e Sequenciamento de Lotes com Linhas de Produção em Paralelo (GLSPPL), proposto por Meyr (2002) e Meyr & Mann (2013), para apresentar três modelos que tratam do planejamento de colheita de safra de cana-de-açúcar considerando as frentes de colheita. De acordo com a analogia utilizada pelos autores, os blocos de colheita são representados pelos lotes do GLSPPL, enquanto as frentes de colheita são representadas pelas linhas de produção.

Diferentemente dos demais trabalhos, que utilizaram predominantemente estimativas de teor de sacarose para a avaliação da maturidade da cana-de-açúcar plantada, que orienta o período em que a área pode ser colhida, Junqueira & Morabito (2017) consideraram o período útil de industrialização (PUI) de cada bloco de colheita para determinar a janela de tempo em que determinada área pode ser colhida. Embora sejam indicadores interessantes da receita esperada da matéria-prima, as estimativas de teor de sacarose podem não se repetir na safra subsequente devido a variações climáticas, de adubação do solo, atuação de pragas, doenças e demais fatores biológicos que possam modificar o comportamento da planta. Já o PUI é proveniente de uma avaliação agronômica abrangente por parte da equipe técnica da usina, bem como das instituições desenvolvedoras de variedades de cana-de-açúcar. Além da maturação, o PUI é formado considerando também outras características da variedade plantada, como a facilidade de brotação, a tendência de florescimento, a resposta a maturadores e à irrigação.

Os modelos apresentados em Junqueira & Morabito (2017) também consideram o potencial de colheita e transporte variável por bloco de colheita. Características como a produtividade esperada do canavial (em toneladas por hectare), o espaçamento de plantio e a quantidade de manobras variam, dependendo da área, e influenciam diretamente na taxa de produção horária de colhedoras. Similarmente, a distância das áreas de colheita à usina, bem como a velocidade de trajeto dos veículos impactam na taxa de transporte dos caminhões. A busca do equilíbrio das capacidades de colheita e transporte é fundamental para a redução de custos com a ociosidade de recursos que, de acordo com os autores, são significativos. Entretanto, a redução de recursos não deve implicar em falha no suprimento de matéria-prima para a indústria, nem em colheita fora da época prevista no PUI, pois esses custos marginais podem ser significativamente mais altos que os de ociosidade. Os demais modelos reportados na literatura consideram apenas a distância como variável dependente da área de colheita, que deveria ser equilibrada ao longo da safra.

Apesar dessa abordagem ser promissora, Junqueira & Morabito (2017) resolveram o problema apenas para exemplares de pequeno porte, com poucas centenas de restrições e variáveis, com dezenas dessas variáveis discretas. Para tratar o planejamento de colheita de safra de cana-de-açúcar de porte real, Higgins et al. (2004b), Higgins & Laredo (2006) bem como Jena & Poggi (2013) propuseram métodos de agregação de blocos de colheita para reduzir a quantidade de variáveis envolvidas. Paralelamente, métodos heurísticos baseados em PIM, tais como relax-and-fix, apresentados por Pochet & Wolsey (2006), e fix-and-optimize, propostos por Helber & Sahling (2010), vêm sendo largamente utilizados para resolver problemas de dimensionamento e sequenciamento de lotes de produção do tipo GLSPPL de grande porte, como, por exemplo, em Beraldi et al. (2008), Ferreira et al. (2009, 2012), Toso et al. (2009) e Helber & Sahling (2010).

Dessa forma, este trabalho utiliza um modelo de PIM para representar o problema de programação e sequenciamento de frentes de colheita de cana-de-açúcar e propõe métodos heurísticos para sua resolução em situações reais, envolvendo três heurísticas integradas: uma heurística de agregação de blocos de colheita, uma heurística construtiva baseada na heurística de PIM relax-and-fix e uma heurística de melhoria baseada na heurística de PIM fix-and-optimize. Os autores não têm conhecimento de outros trabalhos na literatura que apresentaram abordagens de otimização para esse problema nessa linha de pesquisa.

Assim, a seção 2 detalha o PIM para representar esse problema, a seção 3 detalha os métodos heurísticos propostos para a resolução desse problema no planejamento de colheita de safra de cana-de-açúcar, a seção 4 compara o desempenho dos métodos heurísticos propostos, a seção 5 ilustra e analisa em detalhes uma solução do problema e, por fim, a seção 6 apresenta as considerações finais deste estudo e perspectivas para pesquisa futura.

2 Modelo matemático

Na modelagem do problema, considera-se que os blocos de colheita $j=1\mathrm{..}B$, que possuem uma produção estimada de $p_j$, devem ser sequenciados nas frentes de colheita $l=1\mathrm{..}F$ em um horizonte de tempo finito, dividido em macroperíodos $t=1\mathrm{..}P$. Cada macroperíodo pode ainda ser subdividido em microperíodos não sobrepostos, com tamanhos variáveis. Cada frente de colheita possui uma definição específica de subperíodos. O conjunto $S_t$ descreve os s microperíodos em que s = 1,...,${\displaystyle \sum }_{t=1}^T\left|S_t\right|=N$ (vide Figura 1).

Figura 1
Relacionamento entre macro e microperíodos.

O conjunto $B{l}_{jl}$ representa os blocos j que podem ser colhidos pela frente de colheita l. O conjunto $B{s}_{jt}$ representa a janela de tempo, períodos t, em que o bloco j pode ser colhido. O tipo de colheita $m\in M$ pode ser manual ou mecanizada, sendo $F_m$ o conjunto das frentes de colheita que pertencem ao conjunto M. Assim, em cada bloco j, $f_{jm}$ representa a fração da produção do bloco com o tipo de colheita m.

Em todo macroperíodo t, as frentes de colheita devem atender à faixa de demanda da usina $\left[\min d_{mt},\max d_{mt}\right]$ para cada tipo de colheita m. Além disso, há uma quantidade mínima de área a ser liberada, $vi{n}_t$, proveniente do conjunto $V_j$, que representa o conjunto dos blocos que possuem a fração $f{i}_j$ possível de ser irrigada por vinhaça.

Cada frente l possui $N{m}_l$ máquinas (colhedoras ou carregadoras). A capacidade de colheita do bloco j é dada por $co{l}_j$. As máquinas do tipo de colheita m trabalham $H{t}_m$ horas por dia e possuem uma disponibilidade de trabalho idêntica para as l frentes de $K_t$ horas por macroperíodo t. O modal de transporte considerado é totalmente rodoviário, com frota homogênea de $Nt$ de caminhões, que trabalham Htt horas por dia. A frota de veículos possui também uma disponibilidade de trabalho de $K_t$ horas por macroperíodo t. Além disso, um caminhão pode ir para qualquer frente de colheita, seguindo uma regra de despacho dinâmico. Cada bloco j possui uma capacidade de transporte associada $trans{p}_j$, que é função da distância da usina e das condições da malha viária.

A perda de capacidade da frente de colheita durante suas mudanças de bloco é medida por meio do tempo de deslocamento, $s{t}_{ij}$, entre o bloco i, de onde está partindo, ao bloco j, para onde está se deslocando, a qual é calculada a partir da distância entre esses blocos, $dis{t}_{ij}$. A quantidade de veículos com prancha, Np, que possibilitam o transporte das colhedoras (geralmente de esteira) a longas distâncias também influencia na velocidade da mudança de área. Para haver uma mudança dessa não é necessário que toda cana-de-açúcar do bloco seja colhida, entretanto deve-se colher uma quantidade mínima, em toneladas, para o bloco j e a frente l, $b{m}_{lj}$, para que se justifique o deslocamento da frente.

Dentre os modelos apresentados por Junqueira & Morabito (2017), neste estudo considerou-se o chamado Modelo 1B. Esse modelo possui uma solução factível trivial, que é simplesmente não produzir (i.e., não colher nem transportar), pagar pela falha no suprimento e por toda a matéria-prima deixada para ser colhida na safra seguinte, denominada cana bisada, conforme jargão do setor. Com essa característica, soluções factíveis iniciais podem ser facilmente obtidas, potencializando a utilização de métodos heurísticos.

A seguir, esse modelo é resumidamente apresentado, para facilitar a compreensão do leitor dos métodos de solução propostos na próxima seção e baseados nesse modelo, e para que o material do artigo fique autocontido. A Tabela 1 apresenta as variáveis de decisão. Para mais detalhes desses elementos do modelo, o leitor pode consultar Junqueira & Morabito (2017).

Tabela 1
Variáveis de decisão do modelo.

Símbolo	Definição
$x_{ljs}$	produção do bloco j no microperíodo s pela frente l (em toneladas)
$y_{ljs}$	$\{\begin{array}{c}\mathrm{1,}seafrentelestánoblocojnomicroperíodos\\ \mathrm{0,}casocontrário\end{array}$
$z_{lijs}$	$\{\begin{array}{c}\mathrm{1,}seafrentelsedirigedonóiparaojnomicroperíodos\\ \mathrm{0,}casocontrário\end{array}$
$w{m}_{mt}$	toneladas de moagem perdidas do tipo de colheita m no período t
$w{b}_j$	toneladas de matéria-prima deixadas no bloco j

O modelo de PIM é definido por:

$Minmo{\displaystyle \sum }_{t=1}^{T}w{m}_{mt}+bs{\displaystyle \sum }_{j=1}^{B}w{b}_j+md{\displaystyle \sum }_{l=1}^F{\displaystyle \sum }_{j=1}^B{\displaystyle \sum }_{j=1}^{B}dis{t}_{ij}{z}_{lij}$ (1)

Sujeito à:

${\displaystyle \sum }_{l\in F_m}{\displaystyle \sum }_{j=1}^B{\displaystyle \sum }_{s\in S_t}{x}_{ljs}+w{m}_{mt}\ge \min d_{mt}\left(t=\mathrm{1,}\dots ,T\right),\left(m=man,mec\right)$ (2)

${\displaystyle \sum }_{l\in F_m}{\displaystyle \sum }_{j=1}^B{\displaystyle \sum }_{s\in S_t}{x}_{ljs}\le \max d_{mt}\left(t=\mathrm{1,}\dots ,T\right),\left(m=man,mec\right)$ (3)

${\displaystyle \sum }_{l\in F_m}{\displaystyle \sum }_{s=1}^N{x}_{ljs}+w{b}_j=p_j{f}_{jm}\left(j=\mathrm{1,}\dots ,B\right),\left(m=man,mec\right)$ (4)

${\displaystyle \sum }_{m=man}^{mec}{\displaystyle \sum }_{l\in F_m}{\displaystyle \sum }_{j\in V_j}{\displaystyle \sum }_{s\in S_t}f{i}_j\frac{x_{ljs}}{TC{H}_j}\ge vi{n}_t\left(t=\mathrm{1,}\dots ,T\right)$ (5)

${\displaystyle \sum }_{j=1}^B{\displaystyle \sum }_{s\in S_t}\frac{24}{co{l}_{mj}N{m}_{l}H{t}_m}{x}_{ljs}+{\displaystyle \sum }_{i=1}^B{\displaystyle \sum }_{j=1}^B{\displaystyle \sum }_{s\in S_t}\frac{N{m}_l}{Np}s{t}_{ij}{z}_{lijs}\le K_t\left(m=man,mec\right),\forall l\in F_m,\left(t=\mathrm{1,}\dots ,T\right)$ (6)

${\displaystyle \sum }_{m=man}^{mec}{\displaystyle \sum }_{l\in F_m}{\displaystyle \sum }_{j=1}^B{\displaystyle \sum }_{s\in S_t}\frac{24}{trans{p}_{mj}NtHtt}{x}_{ljs}\le K_t\left(t=\mathrm{1,}\dots ,T\right)$ (7)

$x_{ljs}\le \min \left(\frac{co{l}_{mj}N{c}_{l}H{t}_m}{24},\frac{trans{p}_{mj}NtHtt}{24}\right)K_t{y}_{ljs}\left(m=man,mec\right),\forall l\in F_m,\left(j=\mathrm{1,}\dots ,B\right),\left(t=\mathrm{1,}\dots ,T\right),\forall s\in S_t$ (8)

$x_{ljs}\ge b{m}_{lj}\left(y_{ljs}-y_{ljs-1}\right)\left(m=man,mec\right),\forall l\in F_m,\left(j=\mathrm{1,}\dots ,B\right),\left(s=\mathrm{1,}\dots ,N\right)$ (9)

${\displaystyle \sum }_{j\in \left(B{s}_{jt}\cap B{l}_{jl}\right)}{y}_{ljs}=1\left(l=\mathrm{1,}\dots ,F\right),\left(t=\mathrm{1,}\dots ,T\right),\forall s\in S_t$ (10)

$y_{lis-1}={\displaystyle \sum }_{j=1}^B{z}_{lijs}\left(l=\mathrm{1,}\dots ,F\right),\left(i=\mathrm{1,}\dots ,B\right),\left(s=\mathrm{2,}\dots ,N\right)$ (11)

${\displaystyle \sum }_{i=1}^B{z}_{lijs}=y_{ljs}\left(l=\mathrm{1,}\dots ,F\right),\left(j=\mathrm{1,}\dots ,B\right),\left(s=\mathrm{1,}\dots ,N\right)$ (12)

$y_{ljs-1}\ge y_{ljs}\left(l=\mathrm{1,}\dots ,F\right),\left(j=\mathrm{1,}\dots ,B\right),\left(t=\mathrm{1,}\dots ,T\right),\forall s\in \left(S_t\setminus S{0}_t\right)$ (13)

$x_{ljs}\ge 0\left(l=\mathrm{1,}\dots ,F\right),\left(j=\mathrm{1,}\dots ,B\right),\left(s=\mathrm{1,}\dots ,N\right)$ (14)

$y_{ljs}\in \left\{\mathrm{0,1}\right\}\left(l=\mathrm{1,}\dots ,F\right),\left(j=\mathrm{1,}\dots ,B\right),\left(s=\mathrm{1,}\dots ,N\right)$ (15)

$z_{lijs}\ge 0\left(l=\mathrm{1,}\dots ,F\right),\left(i=\mathrm{1,}\dots ,B\right),\left(j=\mathrm{1,}\dots ,B\right),\left(s=\mathrm{1,}\dots ,N\right)$ (16)

$$w{m}_{mt}\ge 0\left(t=\mathrm{1,}\dots ,T\right),\left(m=man,mec\right)$$ (17)

$w{b}_j\ge 0\left(j=\mathrm{1,}\dots ,B\right)$ (18)

A Função Objetivo 1 avalia o impacto dos custos de falha na moagem, cana bisada e deslocamento das frentes de colheita, que também é considerado, embora seja de segunda ordem se comparado aos outros dois. Variáveis de folga medem as toneladas de moagem perdidas, Restrições 2, e para as toneladas de matéria-prima deixadas no campo sem colher, Restrições 4. Para cada um desses parâmetros há um custo por hora ou por tonelada, $mo$ e $bs$, respectivamente. O parâmetro md representa o custo de transporte de uma frente por um quilômetro.

As Restrições 2 a 4 garantem o balanço de massa, ajustando a produção no campo com a demanda da usina. As Restrições 2 e 3 estão relacionadas à demanda de moagem mínima e máxima da usina por período, garantindo que a quantidade de matéria-prima colhida e transportada de um determinado tipo de colheita seja maior que o valor demandado no período. Já as Restrições 4 limitam a colheita e o transporte à disponibilidade de matéria-prima no bloco, estando, assim, relacionadas à produção no campo, respeitando-se o tipo de colheita da frente e do bloco. As Restrições 5 determinam a quantidade mínima de área de vinhaça que deve ser liberada em cada macroperíodo t, respeitando-se a quantidade de área irrigável no bloco, sendo que $TC{H}_j$ representa a produtividade em toneladas por hectare do bloco j.

As Restrições 6 e 7 consideram a capacidade dos recursos de colheita e transporte. As Restrições 6 relacionam o tempo consumido na produção dos recursos de colheita e o tempo consumido no transporte de equipamentos durante a mudança de área com o tempo total disponível no período. O tempo de produção dos recursos leva em conta as características de rendimento do bloco, bem como a adequação agronômica para a colheita no período. O tempo de transporte dos equipamentos na mudança de área leva em conta o tempo gasto nessa mudança, bem como o número de veículos prancha disponíveis para essa operação. As Restrições 7 tratam os recursos de transporte, cujo potencial de produção é considerado por bloco. Nesse caso, não é inserido o tempo gasto com a mudança da frente, pois no momento da mudança de área esses recursos podem transportar a produção das áreas em operação, ou ficam ociosos, consumindo estoque sobre rodas na unidade fabril. Além disso, como não há caminhões específicos por frente de colheita, pois eles podem ir qualquer frente, esse balanço de capacidade pode ser feito apenas para todos os períodos t e não todas as frentes l e períodos t, como na Restrição 6, em que a quantidade de colhedoras é definida por frente de colheita.

As Restrições 8 a 9 acoplam as variáveis $x_{ljs}$ e $y_{ljs}$. As Restrições 8 garantem que quando há produção no bloco j pela frente de colheita l no microperíodo s, a frente está posicionada no mesmo local e momento. Da mesma forma, para a situação oposta, quando a frente não está posicionada, não pode haver produção do bloco. Cabe ressaltar que o limitante superior de $x_{ljs}$ foi considerado o mínimo entre a capacidade de colheita e transporte. Espera-se que a capacidade de colheita seja a mais restritiva, a não ser que uma frente tenha maior capacidade que toda a frota. As Inequações 9 são as restrições de lote mínimo e definem uma quantidade mínima de matéria-prima a ser colhida. Já as Equações 10 garantem que a frente l estará em apenas um bloco j no microperíodo s. Por meio do parâmetro $B{s}_{jt}$, essas restrições impõem que ocorra a colheita nos blocos j permitida pelas janelas de tempo, ou seja, nos macroperíodos t. Além disso, por meio do parâmetro $B{l}_{jl}$, essas equações definem os blocos j que uma frente de colheita l pode colher. Isso ocorre caso seja necessária a setorização de alguma frente. Cabe ressaltar que são previamente fixados $y_{ljs}=0$ para os microperíodos s que estão fora da janela de tempo do bloco j ($B{s}_{jt}$), bem como para as frentes l que não estão aptas a colher ($B{l}_{jl}$).

As Restrições 11 e 12 definem o movimento da frente por meio da variável $z_{lijs}$, consistente com as posições da frente no microperíodo s ($y_{ljs}$) e no microperíodo anterior s-1 ($y_{lis-1}$). As Restrições 13, embora sejam redundantes para a obtenção de uma solução ótima factível, estabelecem que os microperíodos ociosos ocorram apenas no final de cada macroperíodo e, com isso, eliminam soluções equivalentes (simetrias). As Restrições 14, 15, 16, 17 e 18 definem as variáveis não negativas $x_{ljs}$, $z_{lijs}$, $w{m}_{mt}$ e .${{\displaystyle w}}_{bj}$., bem como a variável binária $y_{ljs}$.

3 Métodos de solução propostos

A resolução exata de um problema de porte real com base no modelo utilizado neste estudo pode não ser viável devido ao tamanho do problema. No caso da empresa estudada, que é uma usina de médio porte, ela possuía para a safra estudada cinco frentes de colheita, oito meses de safra (macroperíodos), com 200 dias efetivos de safra e um tempo de permanência médio de três dias da frente por bloco. Portanto, 80 microperíodos e 330 blocos de colheita, que resultariam em 132 mil variáveis inteiras.

A primeira etapa do método proposto foi de agregação dos blocos de colheita, de tal forma que não se descaracterizassem variáveis como a distância do bloco em relação à usina e o potencial de colheita. A segunda etapa envolveu proposição de métodos heurísticos baseados em PIM, como relax-and-fix, para construir uma boa solução inicial factível. A terceira etapa envolveu a utilização do método fix-and-optimize, para buscar a melhoria da solução construída na etapa 2. A Figura 2 ilustra o método heurístico proposto.

Figura 2
Método heurístico proposto.

Uma simples heurística alternativa seria resolver o problema original utilizando um pacote de otimização para resolver o modelo dentro de um tempo máximo limitado, mesmo sem chegar à solução ótima do problema. Nesse caso, também poderia ser adequado utilizar a heurística agregação de blocos (etapa 1). Outra opção seria resolver apenas as etapas 1 e 2.

3.1 Agregação de blocos

Primeiramente, o bloco agregado não pode desconsiderar as janelas de tempo de cada bloco de colheita, representadas por $B{s}_{j\text{'}t}$ no modelo. Essas janelas também podem ser caracterizadas por permitir (1) ou não (0) a colheita em determinado mês (macroperíodo). Se a safra for dividida em oito macroperíodos, por exemplo, cada bloco poderia ter um padrão temporal do tipo |0|0|0|1|1|1|0|0|. Esse código representa que o bloco não pode ser colhido nem nos três primeiros macroperíodos, nem nos dois últimos. A colheita é permitida do quarto ao sexto macroperíodo. Assim, pode haver P padrões temporais, cada um deles caracterizado por permitir (1) ou não (0) a colheita em determinado macroperíodo t ($P{T}_{pt}$). Dessa forma, $B{T}_p$ representa o conjunto dos blocos j’ pertencentes ao padrão temporal p em que $B{s}_{j\text{'}t}$=$P{T}_{pt}$.

Neste estudo, 10 padrões foram detectados, conforme mostra a Tabela 2. O padrão 10 refere-se às áreas ainda a serem plantadas. Os padrões 1 e 2 são de áreas plantadas com variedades tardias, cuja recomendação de colheita é nos últimos dois meses de safra. Os padrões 3, 4, 5 e 6 são de áreas plantadas com variedades médias e os padrões 7, 8 e 9 são de áreas plantadas com variedades precoces, cuja recomendação é colher nos três primeiros meses de safra. Cabe ressaltar que o padrão |1|1|1|1|1|1|1|1| é aquele que possui maior grau de liberdade para encaixar a programação da área em qualquer período, podendo, assim, atuar como fator de equilíbrio das capacidades de colheita e transporte. Após determinado o período ideal de colheita para a área, a equipe agrícola buscar plantar variedades que se adequem ao ambiente de produção e cujo PUI seja compatível.

Tabela 2
Padrões de período (I A1).

Padrão		P1	P2	P3	P4	P5	P6	P7	P8
1	|0|0|0|0|0|0|1|1|	0	0	0	0	0	0	1	1
2	|0|0|0|0|1|1|1|1|	0	0	0	0	1	1	1	1
3	|0|0|0|1|1|0|0|0|	0	0	0	1	1	0	0	0
4	|0|0|0|1|1|1|0|0|	0	0	0	1	1	1	0	0
5	|0|0|1|1|1|1|0|0|	0	0	1	1	1	1	0	0
6	|0|0|1|1|1|1|1|0|	0	0	1	1	1	1	1	0
7	|0|1|1|1|0|0|0|0|	0	1	1	1	0	0	0	0
8	|1|1|1|0|0|0|0|0|	1	1	1	0	0	0	0	0
9	|1|1|1|1|0|0|0|0|	1	1	1	1	0	0	0	0
10	|1|1|1|1|1|1|1|1|	1	1	1	1	1	1	1	1

Além do padrão temporal, os blocos de colheita usualmente considerados nas usinas não ultrapassam em termos espaciais os limites da fazenda, porém, se são conhecidas as coordenadas geográficas do bloco j’ ($X_{j\text{'}},Y_{j\text{'}}$), é possível extrapolar esses limites e tratar uma região geográfica mais ampla. No caso da empresa estudada, a área operacional da usina foi separada em reticulados de 10 km, resultando em um envoltório de 68 km de largura e 71 km de altura. Assim, a quadrícula q é definida pelo conjunto de coordenadas {($X_{q\min },Y_{q\min }$), ($X_{q\min },Y_{q\max }$), ($X_{q\max },Y_{q\min }$), ($X_{q\max },Y_{q\max }$)}, dado que o ponto ($X_{q\min },Y_{q\min }$) representa a coordenada mínima e o ponto ($X_{q\max },Y_{q\max }$) representa a coordenada máxima dessa quadrícula. Para que o bloco j’ pertença à quadrícula q, as seguintes condições devem ser satisfeitas: $X_{q\min }\le X_{j\text{'}}<X_{q\max }$ e $Y_{q\min }\le Y_{j\text{'}}<Y_{q\max }$. Assim, o conjunto $B{E}_q$ agrupa os blocos j’ pertencentes à quadrícula q.

O Gráfico 1 mostra a distribuição geográfica dos blocos de colheita originais e dos agregados. A usina está localizada na coordenada (0,0). Os blocos originais estão em cor cinza. No caso dos blocos agregados, cada tipo de marcador representa um padrão temporal. As cores dos blocos agregados distinguem os diferentes padrões pertencentes a um mesmo tipo de variedade.

Gráfico 1
Distribuição geográfica dos blocos de colheita.

Analisando-se o Gráfico 1, percebe-se que embora haja alguns pontos em que o bloco original coincida com o bloco agregado (quadrícula X=[-40,-30] e Y=[0,10]), especialmente próximo à usina, percebe-se que vários blocos foram agrupados em um (quadrícula X=[-10,0] e Y=[0,10]). Também nessa última quadrícula observa-se que há vários padrões de período, contudo eles não se repetem dentro dela. Há alguns outros casos em que há vários padrões dentro de uma mesma quadrícula (por exemplo, X=[10,20] e Y=[-20,-10]).

Assim, o bloco agregado j, então, é formado pelo conjunto de blocos j’ em que $X_{q\min }\le X_{j\text{'}}<X_{q\max }$, $Y_{q\min }\le Y_{j\text{'}}<Y_{q\max }$ e $B{s}_{j\text{'}t}$=$P{T}_{pt}$, ou seja, $B{A}_j=B{A}_{qp}=B{E}_q\cap B{T}_p$. Por exemplo, para j’ = 1,...,10, q = $\alpha$, $\beta$ e p = A, B, tem-se $B{E}_{\alpha }=\left\{\mathrm{1,3,7,9}\right\}$, $B{E}_{\beta }=\left\{\mathrm{2,4,5,6,8,10}\right\}$, $B{T}_A=\left\{\mathrm{2,4,7,8,9,10}\right\}$, $B{T}_B=\left\{\mathrm{1,3,5,6}\right\}$. Com isso, são formados quatro blocos, j = 1,2,3,4, ou seja, $B{A}_1=B{E}_{\alpha }\cap B{T}_A=\left\{\mathrm{7,9}\right\}$, $B{A}_2=B{E}_{\alpha }\cap B{T}_B=\left\{\mathrm{1,3}\right\}$, $B{A}_3=B{E}_{\beta }\cap B{T}_A=\left\{\mathrm{2,4,8,10}\right\}$, $B{A}_4=B{E}_{\beta }\cap B{T}_B=\left\{\mathrm{5,6}\right\}$.

Uma vez formado o bloco agregado, é necessário determinar suas características agregadas baseadas na unidade de área j’, que seriam: a produção, o potencial de colheita, de transporte e as coordenadas geográficas. A produção do bloco agregado seria a soma de $p_{j\text{'}}$, bem como $co{l}_{mj\text{'}}$, $trans{p}_{mj\text{'}}$, $X_{j\text{'}}$ e $Y_{j\text{'}}$ são determinados por meio da média ponderada pela produção, conforme as Equações 19 a 23.

$p_j={\displaystyle \sum }_{j\text{'}\in B{A}_j}{p}_{j\text{'}}$ (19)

$co{l}_{mj}=\frac{{{\displaystyle \sum }}_{m=man}^{mec}{{\displaystyle \sum }}_{j\text{'}\in B{A}_j}co{l}_{mj\text{'}}*p_{j\text{'}}}{{{\displaystyle \sum }}_{j\text{'}\in B{A}_j}{p}_{j\text{'}}}$ (20)

$trans{p}_{mj}=\frac{{{\displaystyle \sum }}_{m=man}^{mec}{{\displaystyle \sum }}_{j\text{'}\in B{A}_j}trans{p}_{mj\text{'}}*p_{j\text{'}}}{{{\displaystyle \sum }}_{j\text{'}\in B{A}_j}{p}_{j\text{'}}}$ (21)

$X_j=\frac{{{\displaystyle \sum }}_{j\text{'}\in B{A}_j}{X}_{j\text{'}}*p_{j\text{'}}}{{{\displaystyle \sum }}_{j\text{'}\in B{A}_j}{p}_{j\text{'}}}$ (22)

$Y_j=\frac{{{\displaystyle \sum }}_{j\text{'}\in B{A}_j}{Y}_{j\text{'}}*p_{j\text{'}}}{{{\displaystyle \sum }}_{j\text{'}\in B{A}_j}{p}_{j\text{'}}}$ (23)

No exemplo apresentado no início desta seção foi possível reduzir o número de blocos de 330 para 93 para o caso da empresa estudada. Com isso, o número de variáveis inteiras passa a ser 37.200, ao invés de 132 mil, como no problema original. Entretanto, deve-se ter cuidado ao definir o tamanho da quadrícula, pois, de um lado, quanto maior seu tamanho, maior a redução do número de variáveis do problema, mas, de outro, quadrículas muito grandes podem encobrir o efeito da mudança de área. Nesse sentido, é importante que um planejador que conheça bem as áreas da usina possa auxiliar na determinação dos blocos agregados de colheita.

3.2 Heurística de construção relax-and-fix

Como a heurística relax-and-fix possibilita a geração de uma solução inicial factível (heurística construtiva), sua execução sucede a formação dos blocos agregados. A literatura sugere várias estratégias de decomposição das variáveis de decisão para o método relax-and-fix para a resolução do GLSPPL, dentre elas as que utilizam o período (macroperíodos) e os recursos (frentes de colheita). Também é possível uma composição das duas. Ao se utilizar uma estratégia relax-and-fix com base no período, pode-se pensar em pelo menos duas formas: forward, do primeiro ao último período, ou backward, do último ao primeiro.

Seja tt um índice auxiliar idêntico a t, que representa os macroperíodos. Assim, a Figura 3 mostra um algoritmo relax-and-fix proposto para a resolução do PIM proposto com uma estratégia de decomposição da variável $y_{ljs}$ temporal forward e a Figura 4, o mesmo tipo de algoritmo, com uma estratégia temporal backward.

Figura 3
Algoritmo relax-and-fix temporal forward.

Figura 4
Algoritmo relax-and-fix temporal backward.

Nesses dois casos foi possível reduzir o número de microperíodos dentro de cada subproblema. No exemplo tratado na seção 3.1, os 80 microperíodos são reduzidos para 10 em cada subproblema. Com isso, a quantidade de variáveis inteiras passa a ser 4.650, aproximadamente 88% (37.200 variáveis) menor do que o problema com blocos agregados e 97% menor do que o problema original (132 mil variáveis). Apesar de ser necessário resolver o subproblema oito vezes, esse porte para o GLSPPL é significativamente mais tratável computacionalmente.

Seja ll um índice auxiliar idêntico a l, que representa as frentes de colheita (recursos). Assim, a Figura 5 mostra um algoritmo relax-and-fix proposto para a resolução do modelo de programação matemática apresentado com uma estratégia de decomposição da variável $y_{ljs}$ por frente de colheita. No exemplo tratado na seção 3.1, as cinco frentes são reduzidas para uma em cada subproblema. Com isso, a quantidade de variáveis inteiras passa a ser 7.440, aproximadamente 67% (37.200 variáveis) menor do que o problema com blocos agregados e 95% menor do que o problema original (132 mil variáveis). Apesar de ser necessário resolver o subproblema cinco vezes, esse porte para o GLSPPL também é significativamente mais tratável computacionalmente. Comparando-se com a estratégia anterior, o porte do subproblema para a partição temporal é menor, mas a partição por frentes permite que o problema seja resolvido um menor número de vezes.

Figura 5
Algoritmo relax-and-fix por frente de colheita.

3.3 Heurística de melhoria fix-and-optimize

Após a execução do método relax-and-fix, espera-se que seja encontrada uma boa solução inicial factível. Para possibilitar sua melhoria, propõe-se a utilização do método fix-and-optimize. Assim como no relax-and-fix, o fix-and-optimize utiliza uma decomposição da variável $y_{ljs}$ para resolver o GLSPPL com base nas estratégias descritas anteriormente: temporal, por frente (recurso) ou por bloco (produto). Em especial, no fix-and-optimize, as técnicas de overlapping são interessantes, pois as variáveis que não forem inteiras e fizerem parte do subproblema estarão com valores fixos. Com a sobreposição, é possível procurar por soluções diferentes das que foram obtidas na heurística construtiva. A Figura 6 mostra um algoritmo para o método de melhoria fix-and-optimize temporal forward com overlapping de dois períodos. No exemplo tratado na seção 3.1, o porte desse subproblema seria semelhante ao dobro do caso do relax-and-fix com partição temporal das variáveis inteiras, ou seja, teria 9.300 variáveis e deveria ser resolvido sete vezes.

Figura 6
Algoritmo fix-and-optimize temporal forward com overlapping.

4 Comparação dos métodos utilizados

Os métodos heurísticos propostos nas seções 3.1, 3.2 e 3.3, desenvolvidos para o modelo da seção 2, foram combinados e testados nesta seção utilizando-se dados do planejamento de safra da empresa estudada. Assim, a seção 4.1 apresenta a empresa estudada, detalhando os dados de entrada dessa empresa, e a seção 4.2 apresenta os experimentos realizados no estudo.

4.1 Empresa estudada e dados coletados

A usina estudada está localizada no oeste do Estado de São Paulo, Brasil, última área de expansão no estado, e sua fundação aconteceu há menos de 10 anos. Possuiu para a safra 2013/2014 uma moagem de 10.500 t/dia, 2.100.000 t na safra. Em termos de estrutura fundiária, é formada por grandes fazendas, possibilitando a operação com com grandes blocos de colheita não muito dispersos geograficamente. A topografia é caracterizada por suaves ondulações, proporcionando elevado nível de mecanização da área. Assim, a colheita é totalmente mecanizada e a estrutura de colheita também é toda própria.

Dados foram coletados nessa usina num período de aproximadamente três anos, por meio de um trabalho de consultoria com foco na logística de colheita e transporte de um dos autores. Durante esse período, foram obtidas informações criteriosas, na época previstas para a então futura safra de 2013/2014, sobre o potencial de colheita e transporte, capacidade de moagem, quantidade de matéria-prima por área, disponibilidade de tempo dos recursos e quantidade de recursos. Esse levantamento foi realizado em meados de 2012.

Nesta seção são apresentados os principais aspectos referentes aos dados de entrada do modelo aplicado à empresa estudada após a aplicação das heurísticas de agregação (seção 3.1), denominado exemplar A1 – maiores detalhes podem ser vistos em Junqueira (2014). Primeiramente, que todos os blocos de colheita são mecanizáveis e que a colheita será feita com três frentes de cinco colhedoras e duas frentes de quatro. Descontadas as horas paradas para manutenção, abastecimento, troca de turno, horário de refeição e outras paradas operacionais, as colhedoras devem trabalhar 14 horas por dia, segundo a meta adotada pela empresa para a safra. O mesmo vale para as horas trabalhadas dos caminhões, que foi definida pela empresa como 16,6 horas trabalhadas como meta.

Para ser mais fácil a visualização das características dos blocos, a Tabela 3 mostra as estatísticas descritivas básicas das características dos blocos agregados.

Tabela 3
Estatística descritiva das características dos blocos agregados.

	Tamanho do bloco (t)	Col. (t/h)	Transp. (t/h)
Mínimo	267	7	19
Média	22.492	33	28
Máximo	177.583	60	46

O Gráfico 2 complementa essa caracterização, apresentando um histograma do tempo de permanência da frente nos blocos de colheita. Cabe ressaltar que todos os histogramas deste trabalho contêm um gráfico de barra para cada classe do eixo x com a frequência de ocorrência da variável analisada quantificada no eixo y principal (neste caso, o número de blocos) e um gráfico de linha representando o percentual acumulado das classes quantificado no eixo y secundário. Já o Gráfico 3a mostra o histograma do potencial de colheita e o Gráfico 3b, o potencial de transporte.

Gráfico 2
Histograma do tempo de permanência das frentes.

Gráfico 3
Histograma do potencial de colheita (a) e transporte (b).

Analisando-se o tamanho do bloco na Tabela 3 e o tempo de permanência da frente no Gráfico 2, percebe-se que embora haja alguns blocos pequenos, ou seja, com tempo de permanência inferior a um dia, não são significativos, pois representam aproximadamente 10% dos blocos da usina e menos de 1% da quantidade de matéria-prima disponível. Por outro lado, há uma considerável quantidade de blocos grandes, superiores a cinco dias de permanência nas frentes, pois representam aproximadamente 30% dos blocos e 82% da matéria-prima disponível. O peso desses blocos grandes desloca a média para um valor de aproximadamente 10 mil toneladas acima da mediana. Essa característica reduz a complexidade para a obtenção de uma solução boa para o problema, pois diminui a necessidade de mudança de área das frentes de colheita, podendo elas ficarem longos períodos em blocos grandes, com perdas de capacidade com o deslocamento das frentes desprezíveis.

O Gráfico 2 também mostra que em menos de 10% dos blocos o tempo de permanência é inferior a um dia e que em cerca de 30% esse tempo é inferior a três dias. Sendo assim, a quantidade de microperíodos adotados para esse exemplar foi de 10 por frente em um período mensal. Isso equivale a uma mudança de área da frente de colheita a cada três dias. Essa definição é bastante crítica para definir o porte do modelo e, pelos dados apresentados, é bastante razoável pelo fato de três dias de permanência da frente no bloco estar abaixo da média e da mediana.

Ao se analisar o potencial de colheita na Tabela 3 e no Gráfico 3a, percebe-se uma dispersão grande, com potenciais muito baixos, variando entre 7 e 60 toneladas por hora. Nesse caso, o potencial de colheita baixo pode ser atribuído em especial à seca das safras 2010/2011, 2011/2012 e 2012/2013. Para o caso de uma colhedora trabalhando 14 horas por dia, sua produção varia de 100 a 840 toneladas por dia. Para uma produção de 10.500 toneladas por dia, se todas as colhedoras estiverem em uma condição de produção mínima, são necessárias 105 unidades, todavia na melhor condição seriam necessárias apenas 13. Na condição média são necessárias 22,7 colhedoras, sendo que a usina dispõe de 23 no total de suas frentes. As restrições de capacidade de colheita permitem um equilíbrio entre o potencial das áreas disponíveis, a quantidade de colhedoras disponíveis e o tempo de deslocamento das frentes.

Similarmente à colheita, analisando-se o potencial de transporte na Tabela 3 e no Gráfico 3b, percebe-se uma dispersão menor que o potencial de colheita, variando entre 19 e 46 toneladas por hora. Para o caso de um caminhão trabalhando 16,6 horas por dia, sua produção varia de 318 a 840 toneladas por dia. Para uma produção de 10.500 toneladas por dia, se todos os caminhões estiverem em uma condição de produção mínima, seriam necessárias 33 unidades, todavia na melhor condição são necessárias apenas 13. Na condição média seriam necessários 22,7 caminhões, sendo que a usina dispõe de 25 no total. Também as restrições de capacidade de transporte permitem um equilíbrio entre o potencial das áreas disponíveis e a quantidade de caminhões disponíveis.

As distâncias entre blocos foram estimadas com base na distância Euclidiana entre suas coordenadas. Essas distâncias Euclidianas foram corrigidas por um percentual de 30%. Como esse deslocamento forma uma matriz simétrica, foi considerado apenas um lado da matriz e excluiu-se a diagonal, que possui apenas valores iguais a zero.

A partir das distâncias entre os blocos foi calculado o tempo de deslocamento, estimando-se uma velocidade de 40 quilômetros por hora dos veículos de transporte para ida, um tempo de embarque e desembarque das colhedoras de 30 minutos, bem como uma eficiência de 85% nessa operação. Essa eficiência representa os casos em que há mais colhedoras que veículos prancha e em que é necessário que a colhedora aguarde a viagem de retorno de outra colhedora para ser transportada. Assim, o Gráfico 4 mostra o histograma do tempo de deslocamento entre blocos de colheita.

Gráfico 4
Histograma do tempo de deslocamento entre os blocos.

Por fim, a Tabela 4 mostra os parâmetros analisados por macroperíodo. Para a safra foram definidos oito macroperíodos relacionados ao mês com turnos de 24 horas.

Tabela 4
Definição de parâmetros por macroperíodo.

Macroperíodo	Capacidade (h)	Demanda mínima (t)	Demanda máxima (t)	Demanda esperada (t)
P1	496	212.717	221.399	217.058
P2	632	271.144	282.211	276.678
P3	634	271.659	282.747	277.203
P4	684	293.474	305.452	299.463
P5	684	293.474	305.452	299.463
P6	684	271.659	282.747	277.203
P7	632	271.144	282.211	276.678
P8	384	164.642	171.362	168.002

A capacidade se refere às horas disponíveis para as frentes de colheita e o transporte, levando em conta os dias efetivos de safra, ou seja, descontando-se o tempo esperado de parada da usina para manutenção. Já a demanda esperada foi calculada considerando as horas disponíveis no macroperíodo e a moagem horária, que, no caso de 10.500 toneladas por dia, é de 437,5 toneladas por hora. As demandas mínimas e máximas foram calculadas para o caso do modelo, com uma margem de 2% acima ou abaixo da moagem esperada.

Conhecendo-se as demandas de moagem por macroperíodo e conhecendo-se a quantidade de matéria-prima por padrão temporal é possível fazer um balanceamento prévio para identificar infactibilidade na resolução dos modelos propostos. Na Tabela 5 foi feito esse balanceamento preliminar com uma heurística simples. Primeiro, atenderam-se os primeiros macroperíodos, aqueles com os padrões temporais mais precoces, depois os médios e os tardios. Se em algum período faltou matéria-prima, utilizou-se o padrão referente às áreas a serem plantadas.

Tabela 5
Balanceamento preliminar.

Período	11111111	11100000	11110000	01110000	00111100	00111110	00011000	00011100	00001111	00000011	Total
1		217.058									217.058
2		106.568	155.844	14.266							276.678
3				16.510	260.293						277.203
4					299.463						299.463
5					149.876	5.277	8.170	13.313	122.827		299.463
6	234.370								42.833		277.203
7	256.186									20.492	276.678
8	168.002										168.002
Total	658.558	323.626	155.844	30.776	710.031	5.277	8.170	13.313	165.660	20.492	2.091.747

Analisando-se a Tabela 5, percebeu-se que há viabilidade em termos de balanceamento entre moagem e matéria-prima disponível levando-se em conta os padrões temporais. Além disso, percebe-se que as áreas a serem plantadas devem suprir falta de matéria-prima média e tardia nos macroperíodos 6, 7 e 8.

4.2 Implementação computacional e análise dos experimentos

O modelo da seção 2 foi implementado em GAMS versão 24.1.3, utilizando-se como solver o CPLEX 12.5.1.0 para resolver o exemplo A1 em um computador de alto desempenho, com processador Intel i7-3770 com oito núcleos e 16GB de memória RAM. Convém lembrar que esse exemplo envolve mais de 30 mil variáveis e 3 milhões de restrições.

A Tabela 6 descreve os experimentos realizados com o exemplar apresentado na seção 4.1, ao qual foi aplicada a heurística de agregação apresentada na seção 3.1, utilizando-se o modelo da seção 2 sem e com as heurísticas construtivas/de melhoria apresentadas nas seções 3.2 e 3.3. Essa tabela está organizada em cinco blocos de colunas. O primeiro bloco diz respeito aos resultados obtidos pela execução do modelo: valor da função objetivo, Gap de otimalidade e tempo de execução. O segundo bloco mostra a configuração de parâmetros do GAMS/CPLEX utilizada em cada experimento: se as heurísticas RINS, Feasibility Pump e Local Branching estavam ligadas; se os processos básicos do CPLEX de Presolve e planos de corte estavam ligados; e o tempo limite de execução (reslim) nos experimentos em que não foram executadas as heurísticas. O terceiro bloco descreve o desenho da heurística construtiva do tipo relax-and-fix, por meio dos seguintes itens: se a partição em subproblemas foi temporal, com base no recurso (frentes de colheita) ou uma combinação das duas; no caso de ser temporal, se a estratégia foi forward ou backward; se em qualquer uma das duas partições possíveis houve overlapping; e o tempo limite de execução de cada subproblema. O quarto bloco trata da configuração da heurística de melhoria do tipo fix-and-optimize, identificando as seguintes configurações: se a partição em subproblemas foi feita com base nos blocos ou temporal forward; se em qualquer uma das duas partições possíveis houve overlapping; e o tempo limite de execução de cada subproblema. O quinto bloco trata uma modificação feita no modelo, tratando a variável z como real ou inteira. Outros parâmetros adotados que não variaram ao longo dos experimentos foram: a escolha automática do CPLEX da quantidade de núcleos do processador; o número máximo de iterações (iterlim) de 10 milhões; a reserva de memória (workmem) de 30GB; e o Gap de otimalidade aceitável zero para o término da execução (optcr).

Tabela 6
Experimentos computacionais com o modelo.

Experimentos	Resultados			Construtiva (Relax&Fix)						Melhoria (Fix&Optimize)										Modelo Var z
	Valor da função objetivo (em milhares)	Gap (%)	Tempo de execução (hh:mm)	RINS	Feasibility Pump	Local Branching	Presolve	Cortes	reslim (1000s)	Temporal	Forward	Backward	Frente	Overlapping	Reslim (1000s)	Bloco	Temporal (forward)	Overlapping	Reslim (1000s)
1	38.926	99,98	14:30	X	X	X	X	X	52,2											R
2	138.196	100,00	14:30				X	X	52,2											R
3	91	92,82	14:31	X			X	X	52,2											R
4	55	88,16	50:00	X			X	X	180											R
5	21.170	99,97	05:23	X			X	X		X	X				10,8					R
6	4.212	99,84	24:13	X			X	X		X		X			10,8					R
7	28.134	99,98	08:15	X			X	X		X		X		X	10,8					R
8	21.170	99,97	08:00	X			X	X		X	X	X			10,8					R
9	21.982	99,97	96:37	X			X	X		X	X	X	X	X	10,8					R
10	4.212	99,84	24:30	X			X	X		X		X			10,8					Int
11	28.315	99,98	09:26	X		X	X	X		X		X			10,8					R
12	28.315	99,98	10:02	X	X		X	X		X		X			10,8					R
13	4.165	99,84	24:10	X		X	X	X		X		X			10,8					R
14	6.005	99,89	23:51			X	X	X		X		X			10,8					R
15	3.758	99,83	23:59				X	X		X		X			10,8					R
16	40.959	99,98	09:37				X			X		X			10,8					R
17	sem solução		06:12			X		X		X		X			10,8					R
18	59	88,92	50:54				X	X		X		X		X	10,8		X	X	21,6	R
19	59	88,92	49:15				X	X		X		X			10,8		X	X	21,6	R
20	159	95,88	24:44				X	X		X		X			10,8	X		X	2,7	R

Primeiramente, o exemplar foi resolvido pelo modelo sem as heurísticas construtivas e de melhoria, isto é, apenas com o modelo limitado por tempo. Os experimentos 1, 2 e 3 comparam os parâmetros do CPLEX, mostrando que as heurísticas de feasibility pump enfatizando qualidade da solução e local branching não ajudaram a melhorar a qualidade da solução, como mostra o experimento 1. Entretanto, a heurística de RINS trouxe resultados relevantes para a qualidade da solução, por exemplo, quando comparamos o experimento 2 com o 3, apesar de seu Gap ainda ser alto. Dentre as três soluções, cabe ressaltar que a qualidade do experimento 3 é razoavelmente boa, pois apesar de o valor da função objetivo ser maior que zero, a quantidade de falta de cana é desprezível. O valor da função objetivo está baseado em sobra de cana para o ano seguinte, ou seja, em cana bisada em volume pouco significativo, 18.232 toneladas, que representariam menos de dois dias efetivos de moagem. Por ter sido o mais bem-sucedido, o experimento 4 é analisado com mais detalhes na seção 5.

Para os próximos experimentos, que envolvem a aplicação de combinações das heurísticas propostas nas seções 3.2 e 3.3 ao modelo, optou-se, então, por utilizar a configuração com RINS, pré-processamento e plano de cortes ligados, com a ênfase do feasibility pump em otimalidade e local branching desligado. Cabe ressaltar que as soluções dos experimentos 18, 19 e 20 são comparáveis às dos experimentos 3 e 4, com boa qualidade, pois a solução de todos atendeu à demanda de moagem e bisou pouca cana – para maior detalhamento dos testes realizados, vide Junqueira (2014).

Quanto aos tempos computacionais, o Gráfico 5 mostra a evolução dos valores da função objetivo dos experimentos 3 e 4 em função do tempo de execução. Observando-se esse gráfico, percebe-se que, com aproximadamente 38 horas de execução, o experimento 4 obteve uma solução boa, podendo indicar um possível momento de parada da execução do GAMS/CPLEX. O mesmo não ocorre com o experimento 3, que obtém o valor final da solução muito próximo do tempo de execução limite.

Gráfico 5
Função objetivo × tempo de execução.

Assim, pode-se dizer que os experimentos 3, 4, 18, 19 e 20 apresentaram soluções com resultados aceitáveis na prática. Essas soluções podem ser consideras boas, pois foram penalizadas apenas por parte da matéria-prima que não foi possível colher no ano corrente, nesses casos não houve falta de matéria-prima. Embora o experimento 4 tenha apresentado a solução de melhor qualidade, o experimento 3 já apresentou uma solução razoavelmente boa em um tempo de execução muito menor (vide Gráfico 5). Já as heurísticas propostas, obtiveram bons resultados quando combinadas (construtiva + melhoria), mas não se pode dizer que foram melhores que o CPLEX com tempo de execução equivalente limitado. A qualidade dos resultados dos experimentos 18 e 19 foi equivalente à do experimento 4, com tempos de execução similares.

5 Análise da solução obtida

A solução analisada nesta seção refere-se à encontrada no experimento 4, utilizando o exemplar A1 estudado na seção anterior e o modelo da seção 2, uma vez que essa solução apresentou melhor qualidade dentre as que foram consideradas boas. A Tabela 7 mostra a análise da moagem. Esse plano atende em todos os macroperíodos a restrição de moagem mínima, todavia nos períodos P1, P2, P7 e P8 a colheita é superior à moagem mínima. O plano também não extrapola a quantidade moída máxima, chegando ao seu limite em P7 e P8.

Tabela 7
Análise da moagem.

Macroperíodo	Produção (t)	Moagem mínima (t)	Acima da moagem mínima (t)	Moagem máxima (t)	Abaixo da moagem máxima (t)
P1	220.603	212.717	7.886	221.399	797
P2	276.252	271.144	5.107	282.211	5.960
P3	271.659	271.659	0	282.747	11.088
P4	293.474	293.474	0	305.452	11.979
P5	293.474	293.474	0	305.452	11.979
P6	271.659	271.659	0	282.747	11.088
P7	282.211	271.144	11.067	282.211	0
P8	171.362	164.642	6.720	171.362	0
Total	2.080.692	2.049.912	30.780	2.133.582	52.890

Do total de 2.091.747 toneladas de cana disponível, sobraram 11.055 toneladas de cana bisada (3.450 t do bloco b37 e 7.604 t de b75). Nessa solução ainda havia capacidade industrial para absorver essa moagem, porém não havia capacidade para os recursos de colheita e transporte para tal. A Tabela 8 mostra a análise da capacidade dos recursos de colheita e transporte. Essa tabela mostra a quantidade de horas disponíveis por macroperíodo, o quanto foi utilizado pelos recursos de transporte e pelos de colheita. No caso da colheita, a utilização da capacidade foi desmembrada em dois pontos: operação, representando o tempo em que os recursos estiveram de fato trabalhando, e deslocamento, representando o tempo perdido com o deslocamento de recursos entre as áreas.

Tabela 8
Análise da capacidade.

	P1	P2	P3	P4	P5	P6	P7	P8	TOTAL
Disponível	496	632	634	684	684	634	632	384	4.781
Transporte	414	523	615	659	616	545	632	384	4.389
% folga	17	17	3	4	10	14	0	0	8
Colheita
Total operação	487	610	617	666	670	620	522	308	4.500
Deslocamento	9,0	22,5	16,1	18,3	14,5	14,0	33,8	16,3	145
Total	496	632	634	684	684	634	556	324	4.645
% deslocamento	2	4	3	3	2	2	6	5	3
% folga	0	0	0	0	0	0	12	16	3

Analisando-se a Tabela 8 percebe-se que os macroperíodos P7 e P8 exigiram mais dos recursos de transporte, enquanto os recursos de colheita foram mais exigidos de P1 a P6. No geral, houve folga dos recursos, 8% para os recursos de transporte e 3% para os de colheita. Do tempo total dos recursos de colheita, cerca de 3% do tempo foi destinado ao deslocamento, sendo ele mais intenso em P7 e P8. A Tabela 5 da seção 4 mostrou a necessidade de se plantar áreas com variedades para o final da safra. Assim, como as áreas de plantio possuem uma produtividade maior que as demais, proporcionam maior capacidade de colheita aos recursos, originando essa folga adicional nos últimos macroperíodos. Com isso, possibilitou atingir o nível máximo de moagem, passando o recurso de transporte a ser restritivo. Embora o potencial de transporte das áreas de plantio seja menor que a média, sua diferença não é significativa, menos de uma tonelada por hora por veículo.

O Gráfico 6 mostra a distribuição geográfica das frentes e o Gráfico 7, o sequenciamento da frente 1, por exemplo.

Gráfico 6
Distribuição geográfica dos blocos de colheita por frente.

Gráfico 7
Sequenciamento da frente 1.

Observando-se os gráficos, percebe-se que para o plano possibilitar o equilíbrio entre a capacidade de colheita e transporte, uma frente de colheita pode visitar várias vezes um mesmo bloco em épocas diferentes. Um bloco também pode ser visitado por várias frentes, desde que se respeite a quantidade mínima a ser colhida em ambas as passagens pelo bloco. Esse conceito de sequenciamento das frentes de colheita, portanto, difere bastante da visão de setorização predominante na prática das usinas, em que os blocos são previamente fixados às frentes de colheita. Cabe ressaltar que se houver alguma restrição referente ao transporte de funcionários das frentes de colheita, o modelo pode incorporá-la. Na verdade, a setorização é uma prática que levaria a boas soluções se fosse considerado apenas o trajeto mínimo das frentes, sem se considerar o equilíbrio das capacidades de colheita e transporte.

Cabe ressaltar que não setorizar, do ponto de vista do tomador de decisão, pode tornar seu planejamento caótico, caso não tenha à sua disposição uma ferramenta como a que está sendo proposta neste estudo. Além disso, ainda haverá dificuldade na resolução do problema de programação das frentes para a safra toda, o que é o foco deste trabalho.

Por fim, a solução obtida no experimento 4 aplicado ao exemplar A1 estudado utilizando-se o modelo, demonstra uma boa qualidade para a programação das frentes de colheita com as premissas assumidas na seção 4. Embora não tenha sido obtida uma solução com Gap = 0, ou seja, comprovadamente ótima, não houve perda de moagem e a quantidade de matéria-prima bisada pode ser considerada irrelevante. Os dados mostram a importância de se utilizar as capacidades de moagem, transporte e colheita, incluido-se a perda de tempo da frente por mudança de área, da forma proposta neste estudo. Embora a setorização apoie o tomador de decisão para obter soluções sem ferramentas como a proposta neste trabalho, no estudo do caso desta empresa ela não é necessária para a obtenção de uma boa solução.

6 Considerações finais

O problema de programação de frentes de colheita visando equilibrar as capacidades de colheita e transporte é complexo e de porte desafiador. O modelo 1B proposto por Junqueira & Morabito (2017) lida com essa complexidade relacionando variáveis que não são tipicamente tratadas na literatura, como a definição das frentes de colheita como resposta do modelo, bem como o equilíbrio das capacidades de colheita e transporte dependentes dos blocos de colheita. Esse estudo anterior mostrou grande potencial de ganho econômico em um exemplar de pequeno porte.

No presente estudo, métodos de agregação, assim como heurísticas baseada em PIM bem-sucedidas na literatura, como o relax-and-fix e o fix-and-optimise, foram implementadas para a resolução do problema da empresa estudada, de grande porte. Várias combinações dos métodos heurísticos propostos neste trabalho, bem como diferentes parametrizações do CPLEX foram testadas, revelando-se capazes de obter boas soluções para o problema real. Destaca-se que além do porte do problema ter sido desafiador, a folga de capacidade era bastante apertada, especialmente em termos de colheita, e que foi possível obter soluções em tempos computacionais aceitáveis para o planejamento de safra. Cabe ressaltar que para a aplicação desse tipo de abordagem na prática um levantamento de dados criterioso é fundamental, incluindo-se a análise de cada bloco de colheita em termos de capacidade de colheita e transporte, bem como a definição dos PUIs dos blocos, considerando-se os fatores agronômicos preponderantes para a equipe agrícola da usina.

A partir desses estudos, pelo menos três perspectivas de pesquisas futuras interessantes se apresentam. A primeira estaria relacionada ao desenvolvimento de métodos heurísticos para solução e utilização dos modelos 1 e 1B, propostos por Junqueira & Morabito (2017). Uma segunda linha de pesquisa seria a aplicação desse modelo a outras usinas com características diferentes da usina aqui estudada em termos da distribuição de blocos de colheita, potencial de transporte e colheita. Por fim, a terceira linha estaria relacionada com a utilização da solução obtida pelo Modelo 1B para orientação de decisões de plantio e colheita, analisando-se e comparando-se os resultados obtidos por essa abordagem com a praticada anteriormente.

Agradecimentos

Os autores agradecem aos revisores anônimos pelos seus úteis comentários e sugestões de revisão, e à empresa de consultoria e à usina de açúcar e álcool pela colaboração durante este estudo. Esta pesquisa também contou com o apoio da Fapesp e CNPq.

• Suporte financeiro: CNPq e FAPESP.

Referências

• Ahumada, O., & Villalobos, J. R. (2009). Application of planning models in the agri-food supply chain: A review. European Journal of Operational Research, 196(1), 1-20. http://dx.doi.org/10.1016/j.ejor.2008.02.014
  » http://dx.doi.org/10.1016/j.ejor.2008.02.014
• Beraldi, P., Ghiani, G., Grieco, A., & Guerriero, E. (2008). Rolling-horizon and fix-and-relax heuristics for the parallel machine lot-sizing and scheduling problem with sequence-dependent set-up costs. Computers & Operations Research, 35(11), 3644-3656. http://dx.doi.org/10.1016/j.cor.2007.04.003
  » http://dx.doi.org/10.1016/j.cor.2007.04.003
• Ferreira, D., Clark, A. R., Almada-Lobo, B., & Morabito, R. (2012). Single-stage formulations for synchronized two-stage lot sizing and scheduling in soft drink production. International Journal of Production Economics, 136(2), 255-265. http://dx.doi.org/10.1016/j.ijpe.2011.11.028
  » http://dx.doi.org/10.1016/j.ijpe.2011.11.028
• Ferreira, D., Morabito, R., & Rangel, S. (2009). Solution approaches for the soft drink integrated production lot sizing and scheduling problem. European Journal of Operational Research, 196(2), 697-706. http://dx.doi.org/10.1016/j.ejor.2008.03.035
  » http://dx.doi.org/10.1016/j.ejor.2008.03.035
• Grunow, M., Günther, H. O., & Westinner, R. (2007). Supply optimization for the production of raw sugar. International Journal of Production Economics, 110(1-2), 224-239. http://dx.doi.org/10.1016/j.ijpe.2007.02.019
  » http://dx.doi.org/10.1016/j.ijpe.2007.02.019
• Helber, S., & Sahling, F. (2010). A fix-and-optimize approach for the multi-level capacitated lot sizing problem. International Journal of Production Economics, 123(2), 247-256. http://dx.doi.org/10.1016/j.ijpe.2009.08.022
  » http://dx.doi.org/10.1016/j.ijpe.2009.08.022
• Higgins, A. (1999). Optimizing cane supply decisions within a sugar mill region. Journal of Scheduling, 2(5), 329-344. http://dx.doi.org/10.1002/(SICI)1099-1425(199909/10)2:5<229::AID-JOS29>3.0.CO;2-L
  » http://dx.doi.org/10.1002/(SICI)1099-1425(199909/10)2:5<229::AID-JOS29>3.0.CO;2-L
• Higgins, A. (2002). Australian sugar mills optimize harvester rosters to improve production. Interfaces, 32(3), 15-25. http://dx.doi.org/10.1287/inte.32.3.15.41
  » http://dx.doi.org/10.1287/inte.32.3.15.41
• Higgins, A., Antony, G., Sandell, G., Davies, I., Prestwidge, D., & Andrew, B. (2004a). A framework for integrating a complex harvesting and transport system for sugar production. Agricultural Systems, 82(2), 99-115. http://dx.doi.org/10.1016/j.agsy.2003.12.004
  » http://dx.doi.org/10.1016/j.agsy.2003.12.004
• Higgins, A., Haynes, M., Muchow, R., & Prestwidge, D. (2004b). Developing and implementing optimised sugarcane harvest schedules through participatory research. Australian Journal of Agricultural Research, 55(3), 297-306. http://dx.doi.org/10.1071/AR03172
  » http://dx.doi.org/10.1071/AR03172
• Higgins, A., & Laredo, L. A. (2006). Improving harvesting and transport planning within a sugar value chain. The Journal of the Operational Research Society, 57(4), 367-376. http://dx.doi.org/10.1057/palgrave.jors.2602024
  » http://dx.doi.org/10.1057/palgrave.jors.2602024
• Higgins, A., & Muchow, R. (2003). Assessing the potential benefits of alternative cane supply arrangements in the Australian sugar industry. Agricultural Systems, 76(2), 623-638. http://dx.doi.org/10.1016/S0308-521X(02)00031-8
  » http://dx.doi.org/10.1016/S0308-521X(02)00031-8
• Higgins, A., Muchow, R., Rudd, A., & Ford, A. (1998). Optimizing harvest date in sugar production: a case study for the Mossman mill region in Australia I. Development of operations research model and solution. Field Crops Research, 57(2), 153-162. http://dx.doi.org/10.1016/S0378-4290(97)00116-0
  » http://dx.doi.org/10.1016/S0378-4290(97)00116-0
• Jena, S. D., & Poggi, M. (2013). Harvest planning in the Brazilian sugar cane industry via mixed integer programming. European Journal of Operational Research, 230(2), 374-384. http://dx.doi.org/10.1016/j.ejor.2013.04.011
  » http://dx.doi.org/10.1016/j.ejor.2013.04.011
• Junqueira, R. A. R. J. (2014). Programação das frentes de colheita de cana-de-açúcar: uma modelagem visando o equilíbrio das capacidades de colheita e transporte (Tese de doutorado). Departamento de Engenharia de Produção, Universidade Federal de São Carlos, São Carlos.
• Junqueira, R. A. R., & Morabito, R. (2012). Production and logistics planning considering circulation taxes in a multi-plant seed corn company. Computers and Electronics in Agriculture, 84, 100-110. http://dx.doi.org/10.1016/j.compag.2012.02.019
  » http://dx.doi.org/10.1016/j.compag.2012.02.019
• Junqueira, R. A. R., & Morabito, R. (2017). Abordagens de otimização para a programação e sequenciamento das frentes de colheita de cana-de-açúcar. Gestão & Produção, 24(2), 407-422.
• Meyr, H. (2002). Simultaneous lotsizing and scheduling on parallel machines. European Journal of Operational Research, 139(2), 277-292. http://dx.doi.org/10.1016/S0377-2217(01)00373-3
  » http://dx.doi.org/10.1016/S0377-2217(01)00373-3
• Meyr, H., & Mann, M. (2013). A decomposition approach for the general lotsizing and scheduling problem for parallel production lines. European Journal of Operational Research, 229(3), 718-731. http://dx.doi.org/10.1016/j.ejor.2013.03.036
  » http://dx.doi.org/10.1016/j.ejor.2013.03.036
• Plà, L. M., Sandars, D. L., & Higgins, A. J. (2014). A perspective on operational research prospects for agriculture. The Journal of the Operational Research Society, 68(7), 1078-1089. http://dx.doi.org/10.1057/jors.2013.45
  » http://dx.doi.org/10.1057/jors.2013.45
• Plà-Aragones, L. M. (2015). Handbook of operations research in agriculture and the agri-food industry. International Series in Operation Research & Management Science, 224 http://dx.doi.org/10.1007/978-1-4939-2483-7
  » http://dx.doi.org/10.1007/978-1-4939-2483-7
• Pochet, Y., & Wolsey, L. A. (2006). Production planning by mixed integer programming New York: Springer Science+Business Media.
• Toso, E. A. V., Morabito, R., & Clark, A. R. (2009). Lot sizing and sequencing optimisation at an animal-feed plant. Computers & Industrial Engineering, 57(3), 813-821. http://dx.doi.org/10.1016/j.cie.2009.02.011
  » http://dx.doi.org/10.1016/j.cie.2009.02.011

Datas de Publicação

• Publicação nesta coleção
  23 Out 2017
• Data do Fascículo
  Jan-Mar 2018

Histórico

• Recebido
  09 Nov 2015
• Aceito
  22 Abr 2016