Faculteit Economie en Bedrijfskunde Vakgroep Operationeel Onderzoek
Academiejaar 2011–2012
Een vergelijkende studie van exacte en meta-heuristische oplossingsmethodes voor het optimaliseren van de planning van ongerelateerde parallelle machines in de textielindustrie
Louis-Philippe Kerkhove
Promotor: Prof. dr. M. Vanhoucke
Scriptie voorgedragen tot het behalen van de graad van Master in de Toegepaste Economische Wetenschappen: Handelsingenieur
PERMISSION Ondergetekende verklaart dat de inhoud van deze masterproef mag geraadpleegd en/of gereproduceerd worden, mits bronvermelding. Louis-Philippe Kerkhove Gent, Mei 2012
Woord vooraf Dit werkstuk kwam tot stand met de hulp van heel wat personen, graag zou ik alle personen die een bijdrage hebben geleverd aan dit werkstuk willen bedanken, waarvan enkele personen in het bijzonder. Vooreerst mijn oprechte dank aan mijn promotor prof. dr. Mario Vanhoucke, prof. dr. Broos Maenhout en dr. Vincent Van Peteghem en de rest van het academisch personeel op de vakgroep operations waar ik steeds terecht kon met mijn vragen. Daarnaast wil ik graag ook alle mensen van het bedrijf bedanken voor hun enthousiaste medewerking aan dit project. In het bijzonder Stefaan, Bart, Petra en Veronique die bijzonder veel van hun tijd hebben willen investeren in dit project. Ook Frederick Heene zou ik graag willen bedanken voor zijn constructieve feedback tijdens de laatste fase van dit project, evenals Nikolaas Messely en Bernard Dierinck voor hun hulp bij het overlezen van dit werkstuk. Ten slotte nog een speciaal woord van dank aan mijn vriendin, mijn ouders en de ouders van mijn vriendin voor de morele steun gedurende dit project. Louis-Philippe Kerkhove Gent 20 mei 2012
ii
Inhoudsopgave Lijst van tabellen
vii
Lijst van figuren
ix
1 Inleiding
1
I
3
Analyse van de probleemstelling
2 Probleemsituering 2.1 Bedrijfseconomische context . . . . . . . . . . . 2.2 Het productieproces . . . . . . . . . . . . . . . 2.2.1 De stappen in het productieproces . . . 2.2.2 Observaties in de productiehal . . . . . 2.2.3 Bepalen van de scope van het probleem 2.3 Gedetailleerde beschrijving van het probleem . 2.4 Huidige planningsmethode . . . . . . . . . . . . 2.5 Analoge problemen in de literatuur . . . . . . . 2.5.1 Situering van de aard van het probleem 2.5.2 Exacte oplossingsmethodes . . . . . . . 2.5.3 (Meta-)Heuristische oplossingsmethodes 2.5.4 Conclusies uit literatuuronderzoek . . .
. . . . . . . . . . . .
4 4 5 5 7 8 9 10 11 11 13 14 15
3 Doelfunctie en inputparameters 3.1 Frequent gebruikte doelfuncties uit de literatuur . . . . . . . . . . . . . . . . 3.2 Logische Samenstelling van de doelfunctie . . . . . . . . . . . . . . . . . . . . 3.3 Wiskundige formulering van de doelfunctie . . . . . . . . . . . . . . . . . . .
17 17 18 20
4 Data generator 4.1 Input voor data generator . . . . . . . . . . 4.1.1 Namen voor illustratieve doeleinden 4.1.2 Productiespecificaties . . . . . . . . 4.1.3 Categorie¨en van kwaliteiten . . . . . 4.1.4 Instellingen van systemen . . . . . . 4.1.5 Nationaliteiten en ABC scores . . .
22 23 23 23 23 24 24
iii
. . . . . .
. . . . . .
. . . . . . . . . . . .
. . . . . .
. . . . . . . . . . . .
. . . . . .
. . . . . . . . . . . .
. . . . . .
. . . . . . . . . . . .
. . . . . .
. . . . . . . . . . . .
. . . . . .
. . . . . . . . . . . .
. . . . . .
. . . . . . . . . . . .
. . . . . .
. . . . . . . . . . . .
. . . . . .
. . . . . . . . . . . .
. . . . . .
. . . . . . . . . . . .
. . . . . .
. . . . . . . . . . . .
. . . . . .
. . . . . . . . . . . .
. . . . . .
. . . . . . . . . . . .
. . . . . .
. . . . . . . . . . . .
. . . . . .
. . . . . . . . . . . .
. . . . . .
. . . . . . . . . . . .
. . . . . .
. . . . . .
Inhoudsopgave
4.2
II
iv
4.1.6 Transportafstanden . . . . . . . . . . . . . 4.1.7 Ordergroottes . . . . . . . . . . . . . . . . 4.1.8 Machinegroepen . . . . . . . . . . . . . . 4.1.9 Centimeters per omwenteling . . . . . . . Procedure voor datageneratie . . . . . . . . . . . 4.2.1 Informatie over orders . . . . . . . . . . . 4.2.2 Machine informatie . . . . . . . . . . . . . 4.2.3 Parameters met betrekking tot productie 4.2.4 Vaste kosten voor laattijdig leveren . . . . 4.2.5 Variabele kosten voor laattijdig leveren . 4.2.6 Maximale waarden voor indices . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
Ontwerpen van oplossingsmethodes
5 Exact optimalisatiemodel 5.1 Modelformulering . . . . . . 5.1.1 Indices . . . . . . . . 5.1.2 Variabelen . . . . . . 5.1.3 Parameters . . . . . 5.1.4 Doelfunctie . . . . . 5.1.5 Randvoorwaarden . 5.1.6 Beknopte duiding bij 5.2 Resultaten . . . . . . . . . . 5.2.1 Breedte onderzoek . 5.2.2 Diepte onderzoek . .
24 25 25 25 25 26 27 27 29 30 30
32
. . . . . . . . . .
33 33 33 34 35 35 35 36 37 38 38
. . . . . . . . . . . .
42 42 43 43 45 45 45 45 45 46 46 50 56
7 Genetische metaheuristiek 7.1 Werking van de metaheuristiek . . . . . . . . . . . . . . . . . . . . . . . . . . 7.1.1 Representatie van productieschemas . . . . . . . . . . . . . . . . . . .
57 58 58
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vergelijkingen . . . . . . . . . . . . . . . . . . . . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
6 Exacte oplossingsmethode: Illustraties en duiding 6.1 Indices . . . . . . . . . . . . . . . . . . . . . . . . . . 6.2 Beslissingsvariabelen . . . . . . . . . . . . . . . . . . 6.2.1 De positie selectie variabele . . . . . . . . . . 6.2.2 De lateness variabele . . . . . . . . . . . . . . 6.2.3 Tardiness variabele . . . . . . . . . . . . . . . 6.2.4 Tijdstip van opleveren . . . . . . . . . . . . . 6.2.5 Vaste kosten voorgesteld door λ . . . . . . . . 6.3 Parameters . . . . . . . . . . . . . . . . . . . . . . . 6.4 Randvoorwaarden opgelegd aan het model . . . . . . 6.4.1 Afdwingen van een werkbaar schema . . . . . 6.4.2 Tijdstip bepalende randvoorwaarden . . . . . 6.4.3 Eigenschappen van variabelen bepalen . . . .
. . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . . . .
Inhoudsopgave
7.2
7.3
7.4
v
7.1.2 Creatie van de populatie . . . . . . . . . . . . . . . . . 7.1.3 Selectie van ouders . . . . . . . . . . . . . . . . . . . . 7.1.4 Generatie van nakomelingen . . . . . . . . . . . . . . . 7.1.5 Mutaties . . . . . . . . . . . . . . . . . . . . . . . . . . 7.1.6 Generatieovergang . . . . . . . . . . . . . . . . . . . . Kalibratie van parameters . . . . . . . . . . . . . . . . . . . . 7.2.1 Omvang van de populatie . . . . . . . . . . . . . . . . 7.2.2 Pressure% . . . . . . . . . . . . . . . . . . . . . . . . . 7.2.3 Kans op mutaties . . . . . . . . . . . . . . . . . . . . . 7.2.4 Aantal generaties . . . . . . . . . . . . . . . . . . . . . 7.2.5 Aantal nakomelingen . . . . . . . . . . . . . . . . . . . 7.2.6 Toepassing op probleem met realistische grootte . . . 7.2.7 Overzicht van parameters . . . . . . . . . . . . . . . . Uitbreiding van de functionaliteiten . . . . . . . . . . . . . . . 7.3.1 Enkelvoudige heuristieken . . . . . . . . . . . . . . . . 7.3.2 Intelligente iteratie limieten voor het aantal generaties Besluit . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
58 59 59 60 61 62 64 65 66 67 68 68 70 71 71 75 76
8 Simulated Annealing Metaheuristiek 8.1 Werking van simulated annealing metaheuristiek . . . . . . . . . . . . . . . . 8.1.1 Basisprincipes van simulated annealing . . . . . . . . . . . . . . . . . . 8.1.2 Methodes voor lokaal zoeken binnen de simulated annealing heuristiek 8.2 Kalibratie van parameters voor simulated annealing metaheuristiek . . . . . . 8.2.1 Schatting van optimale parameterwaarden . . . . . . . . . . . . . . . . 8.2.2 Starttemperatuur . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8.2.3 Eindtemperatuur . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8.2.4 Alfa . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8.2.5 Maximaal aantal lokale zoekbewegingen . . . . . . . . . . . . . . . . . 8.2.6 Kans op willekeurig zoeken . . . . . . . . . . . . . . . . . . . . . . . . 8.2.7 Kans op wissel . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8.2.8 Overzicht van gekalibreerde parameterwaarden voor simulated annealing 8.3 Besluit . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
77 77 78 79 87 88 89 90 91 92 92 93 94 94
III
96
Resultaten en uitbreidingen
9 Performantie van heuristieken 9.1 Vergelijking met exacte oplossingstechnieken . . . . . . . . . 9.2 Toepassing op realistische probleemgrootte . . . . . . . . . . 9.2.1 Performantie met standaardparameters uit kalibratie 9.2.2 Performantie met hogere toegelaten rekentijd . . . . 9.3 Conclusie . . . . . . . . . . . . . . . . . . . . . . . . . . . .
. . . . .
. . . . .
. . . . .
. . . . .
. . . . .
. . . . .
. . . . .
. . . . .
. . . . .
. . . . .
97 97 98 98 99 102
Inhoudsopgave
vi
10 Conflicten door beperkte beschikbaarheid van personeel 103 10.1 Analyse van conflicten bij geoptimaliseerde schema’s . . . . . . . . . . . . . . 103 10.2 Invloed van het aantal machines op het voorkomen van personeelsconflicten . 105 10.3 Besluit . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 107 11 Slotbeschouwing
108
IV
109
Bijlagen
A Uitbreiding op exacte oplossingsmethode: ten A.1 Indices . . . . . . . . . . . . . . . . . . . . A.2 Variabelen . . . . . . . . . . . . . . . . . . A.2.1 Binaire variabelen . . . . . . . . . A.2.2 Integer variabelen . . . . . . . . . A.2.3 Continue Variabelen . . . . . . . . A.3 Parameters . . . . . . . . . . . . . . . . . A.4 Doelfunctie . . . . . . . . . . . . . . . . . A.5 randvoorwaarden . . . . . . . . . . . . . . A.6 Beknopte duiding bij vergelijkingen . . . .
integratie van personeelsconflic110 . . . . . . . . . . . . . . . . . . . . 111 . . . . . . . . . . . . . . . . . . . . 111 . . . . . . . . . . . . . . . . . . . . 111 . . . . . . . . . . . . . . . . . . . . 111 . . . . . . . . . . . . . . . . . . . . 112 . . . . . . . . . . . . . . . . . . . . 112 . . . . . . . . . . . . . . . . . . . . 113 . . . . . . . . . . . . . . . . . . . . 113 . . . . . . . . . . . . . . . . . . . . 115
B Uitbreiding op exacte oplossingsmethode: Illustraties en duiding B.1 Indices . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . B.2 Beslissingsvariabelen . . . . . . . . . . . . . . . . . . . . . . . . . . . B.2.1 De begin- en eindtijd variabelen . . . . . . . . . . . . . . . . B.2.2 Werkblokken toekennen aan omschakeling . . . . . . . . . . . B.2.3 Exacte duurtijd van omschakeling . . . . . . . . . . . . . . . B.3 Parameters . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . B.4 Randvoorwaarden opgelegd aan het model . . . . . . . . . . . . . . . B.4.1 Allocatie van schaarse grondstoffen . . . . . . . . . . . . . . . B.4.2 Eigenschappen van variabelen bepalen . . . . . . . . . . . . . B.5 Problemen bij implementatie van optimalisatiemodel . . . . . . . . . Bibliografie
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
118 118 119 119 120 120 120 121 121 131 132 133
Lijst van tabellen 4.1
Illustratie berekening omschakeltijden . . . . . . . . . . . . . . . . . . . . . .
30
5.1 5.2 5.3 5.4 5.5 5.6
Indices gebruikt in optimalisatiemodel . . . . . . . . . . . . . . . . . . . . . . Binaire variabelen gebruikt in optimalisatiemodel . . . . . . . . . . . . . . . . Integer variabelen gebruikt in optimalisatiemodel . . . . . . . . . . . . . . . . Continue variabelen gebruikt in optimalisatiemodel . . . . . . . . . . . . . . . Parameters gebruikt in optimalisatiemodel . . . . . . . . . . . . . . . . . . . . Gemiddelde relatieve afwijking van de beste oplossing ten opzichte van de best gekende lower bound, voor verschillende combinaties van machines en jobs. .
33 34 34 34 35
6.1 6.2 6.3 6.4 6.5 6.6 6.7 6.8
Illustratie variabele xjkm . . . . . . . . . . . . Voorbeeld parameterinput voor optimalisatie Voorbeeld matrix van omschakeltijden . . . . Illustratie randvoorwaarde 5.2 . . . . . . . . . Illustratie randvoorwaarde 5.3 . . . . . . . . . Illustratie randvoorwaarde 5.4 . . . . . . . . . Illustratie randvoorwaarde 5.5 . . . . . . . . . Input voor illustratie randvoorwaarde 5.6 . .
44 47 47 48 48 49 51 52
7.1 7.2
Basiswaarden voor parameters gebruikt tijdens kalibratie. . . . . . . . . . . . 63 Samenvatting van de waardes van de parameters zoals bepaald door de kalibratie. 71
8.1 8.2 8.3 8.4 8.5 8.6
Parameters van simulated annealing heuristiek . . . . . . . . . . . . . . . . . Tabel met productietijden voor illustratie lokaal zoeken bij SA . . . . . . . . Tabel met effici¨enties voor illustratie lokaal zoeken bij SA . . . . . . . . . . . Volgorde wissels bij intelligent lokaal zoeken in SA . . . . . . . . . . . . . . . Startwaarden voor kalibratie van simulated annealing heuristiek. . . . . . . . Overzicht van de gekalibreerde waarden voor de simulated annealing heuristiek.
A.1 A.2 A.3 A.4 A.5
Indices gebruikt in optimalisatiemodel . . . . . . . Binaire variabelen gebruikt in optimalisatiemodel . Integer variabelen gebruikt in optimalisatiemodel . Continue variabelen gebruikt in optimalisatiemodel Parameters gebruikt in optimalisatiemodel . . . . .
vii
. . . . . . . .
. . . . . . . .
. . . . . . . .
. . . . . . . .
. . . . .
. . . . . . . .
. . . . .
. . . . . . . .
. . . . .
. . . . . . . .
. . . . .
. . . . . . . .
. . . . .
. . . . . . . .
. . . . .
. . . . . . . .
. . . . .
. . . . . . . .
. . . . .
. . . . . . . .
. . . . .
. . . . . . . .
. . . . .
. . . . . . . .
. . . . .
. . . . . . . .
. . . . .
. . . . . . . .
. . . . .
. . . . . . . .
. . . . .
. . . . . . . .
. . . . .
41
79 83 84 88 89 95 111 111 111 112 112
Lijst van tabellen B.1 B.2 B.3 B.4 B.5 B.6 B.7
Voorbeeld matrix van omschakeltijden . . . . . . . . . . Illustratie randvoorwaarde A.11 t.e.m. A.13 . . . . . . . Een correcte oplossing voor de illustratie van A.11 t.e.m. Illustratie randvoorwaarde A.16 . . . . . . . . . . . . . . Illustratie randvoorwaarde A.18 . . . . . . . . . . . . . . Illustratie randvoorwaarde A.20 en A.21 . . . . . . . . . Illustratie randvoorwaarde A.23 . . . . . . . . . . . . . .
viii . . . . . . . . A.13 . . . . . . . . . . . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
122 122 125 127 128 129 131
Lijst van figuren 2.1 2.2
Schematische voorstelling van het breiproces . . . . . . . . . . . . . . . . . . . Gedetailleerde weergave van schema . . . . . . . . . . . . . . . . . . . . . . .
5 9
3.1
Samenstelling van de vaste en variabele kosten in de doelfunctie . . . . . . . .
19
4.1
Schematische weergave van input en output van de data generator . . . . . .
22
5.1
Weergave van rekentijd nodig voor oplossen van verschillende job en machine combinaties met exact optimalisatiemodel, optimalisatie stopgezet indien geen resultaat na 4000 seconden. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39 Exacte oplossingen gevonden binnen 8000 seconden voor het plannen van 15 jobs op M machines, per machinecategorie werden tien datasets getest. . . . . 40 Probleemgroottes die binnen een kort tijdsbestek exact kunnen worden opgelost. 41
5.2 5.3 7.1 7.2 7.3 7.4 7.5 7.6 7.7 7.8 7.9 7.10 7.11 7.12 7.13
Representatie van genetisch algoritme . . . . . . . . . . . . . . . . . . . . . . Ouders geselecteerd uit populatie met behulp van tournament selectie mechanisme. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Nakomelingen gengenereerd uit ouders met willekeurige crossover. . . . . . . Mutaties bij nakomelingen bevorderen de diversiteit in de populatie. . . . . . Kalibratie van populatieomvang . . . . . . . . . . . . . . . . . . . . . . . . . . Kalibratie van pressure% . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Kalibratie van kans op mutatie bij nakomelingen . . . . . . . . . . . . . . . . Kalibratie van aantal generaties . . . . . . . . . . . . . . . . . . . . . . . . . . Kalibratie van aantal nakomelingen . . . . . . . . . . . . . . . . . . . . . . . . Invloed van het aantal generaties op de rekentijd en de waarde van de doelfunctie bij toepassing op realistische probleemgrootte. . . . . . . . . . . . . . Invloed van het aantal nakomelingen op de rekentijd en de waarde van de doelfunctie bij toepassing op realistische probleemgrootte. . . . . . . . . . . . Gemiddelde afwijking ten opzichte van de optimale oplossing voor de verschillende machine selectie mechanismen voor 40 willekeurige parametersets. . . . Gemiddelde afwijking ten opzichte van de optimale oplossing voor de verschillende volgorde toekennende mechanismen voor 40 willekeurige parametersets.
ix
58 59 60 61 64 65 66 67 68 69 70 73 74
Lijst van figuren 7.14 Gemiddelde cumulatieve afwijking ten opzichte van de optimale oplossing voor combinatie van machine en volgorde heuristieken voor 40 willekeurige parametersets. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8.1 8.2 8.3 8.4 8.5 8.6 8.7 8.8 8.9 9.1 9.2 9.3 9.4
Voorbeeldschema ter illustratie van simulated annealing metaheuristiek . . . Bepalen van de volgorde van de wissels bij intelligent zoeken. . . . . . . . . . Gemiddelde kans op accepteren van nieuwe oplossing met eerste schattingen van parameters. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Invloed van de starttemperatuur op de doelfunctie en de rekentijd. . . . . . . Invloed van de eindtemperatuur op de doelfunctie en de rekentijd. . . . . . . Invloed van alfa op de doelfunctie en de rekentijd. . . . . . . . . . . . . . . . Invloed van van Mmax op doelfunctie en rekentijd. . . . . . . . . . . . . . . . Invloed van kans op willekeurig zoeken op doelfunctie en rekentijd. . . . . . . Invloed van kans op wissel op doelfunctie en rekentijd. . . . . . . . . . . . . . Performantie van exacte en meta heuristische oplossingsmethodes getest op probleem dat 15 jobs inplant op 4 machines. . . . . . . . . . . . . . . . . . . . Performantie van metaheuristieken voor verschillende rekentijden. . . . . . . . Evolutie naar oplossing van simulated annealing heuristiek voor verschillende waarden van α . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Evolutie naar oplossing van genetisch algoritme voor verschillende hoeveelheid nakomelingen per generatie. . . . . . . . . . . . . . . . . . . . . . . . . . . . .
10.1 Kans op tekort aan personeel door gelijktijdige omschakelingen op een willekeurig moment. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . B.1 Illustratie van het gebruik van de begin- en starttijd variabelen . . . . . . . . B.2 Aantal variabelen nodig voor het oplossen van het exact optimalisatiemodel, uitgezet op logaritmische schaal . . . . . . . . . . . . . . . . . . . . . . . . . .
x
75 81 87 89 90 91 92 93 94 95
98 100 101 102
106 119 132
Hoofdstuk 1
Inleiding Het onderwerp van dit werkstuk is de optimalistatie van de planning van het breiproces van een Vlaamse textielfabrikant met een internationaal productieapparaat. Het breiproces zoals dit bij dit bedrijf bestaat wordt in de literatuur omschreven als het parallel machine scheduling (PMS) probleem. Het doel van dit werkstuk is om een conceptuele methode aan te reiken waarmee de planning van dit proces op een betere manier kan tegemoetkomen aan de wensen van het management. Dit betekent dat de economische waarde die gecre¨eerd wordt door het productieproces wordt gemaximaliseerd. De motivatie voor dit onderzoek is terug te vinden in het groot aantal laattijdig geleverde orders. Een meer effectief gebruik van het productieapparaat zou moeten leiden tot het tijdiger opleveren van orders. Dit resulteert in meer tevreden klanten, wat een hogere economische waarde betekent voor zowel de klant als het productiebedrijf. Dit werkstuk bestaat uit vier grote delen. Deel I is descriptief en omvat een uitgebreide analyse van het probleem en de doelstellingen van de planning. Dit gebeurt in de eerste plaats door zowel kwalitatieve als kwantitatieve observaties in de productiehal en gesprekken met een groot aantal werknemers van het bedrijf die met dit proces te maken hebben. Het probleem wordt ook gestiueerd in de literatuur, waarbij de nadruk wordt gelegd op de gebruikte oplossingsmethodes. Daarnaast wordt alle data die nodig is om dit probleem te bevatten ganalyseerd en vervat in een datagenerator. In deel II van dit werkstuk worden verschillende oplossingsmethodes ontwikkeld en gekalibreerd. Het gaat hier meer specifiek over een exacte optimalisatie (MIP) en twee meta1
Hoofdstuk 1. Inleiding
2
heuristische technieken. Deel III evalueert de performantie van deze verschillende oplossingsmethodes. Daarnaast wordt ook kort dieper ingegaan op mogelijke conflicten door bepekte beschikbaarheid van personeel voor het omschakelen van machines. Finaal wordt een algemene conclusie geformuleerd. Deel IV omvat de bijlages en bespreekt een theoretisch optimalisatiemodel dat naast de andere productiebeperkingen (zoals productiemogelijkheden...) ook de beschikbaarheid van personeel opneemt in het model. Dit model was helaas te complex om gebruikt te kunnen worden gegeven de beschikbare middelen.
Deel I
Analyse van de probleemstelling
3
Hoofdstuk 2
Probleemsituering Dit hoofdstuk heeft als doel het planningsprobleem grondig te analyseren en te situeren. Dit gebeurt door eerst de bedrijfseconomische context te bespreken, vervolgens wordt het productieproces aan een grondige analyse onderworpen. In het laatste deel van dit hoofdstuk wordt het afgebakende probleem gesitueerd in de literatuur en worden verschillende oplossingsmethodes uit de literatuur vergeleken.
2.1
Bedrijfseconomische context
De onderzoeksvraag van dit werkstuk komt voort uit de bedrijfseconomische problematiek van deze textielfabrikant. Een groeiende vraag naar producten veroorzaakt een zwaardere belasting van de productiemiddelen. Deze situatie maakt de sector eveneens aantrekkelijker voor concurrenten, die meestal een lage-kostenstrategie volgen. De verhoogde belasting van het productieapparaat ligt aan de oorzaak van de verlenging van de leveringstermijnen. Dit gecombineerd met de verhoogde concurrentie van lage-kosten leveranciers zorgt voor een potentieel gevaarlijke situatie. De verliezen veroorzaakt door laattijdig leveren worden steeds significanter. Het productieapparaat van het bedrijf bevindt zich gedeeltelijk in lageloonlanden, maar tevens voor een belangrijk deel op locaties dicht bij de afnemers. Ondanks deze locatie slaagt het bedrijf er meestal niet in om sneller te leveren dan concurrenten die hun producten over veel grotere afstanden moeten transporteren. Een betere planning zou ervoor kunnen zorgen dat deze locatie in het voordeel van het bedrijf kan worden gebruikt. Indien dit moeilijk haalbaar blijkt gegeven de huidige productiemiddelen kunnen vragen worden gesteld bij het 4
Hoofdstuk 2. Probleemsituering
5
voortbestaan van deze (duurdere) productielocaties. Naast het verhoogde economisch risico veroorzaakt het laattijdig leveren ook verhoogde kosten voor levering en administratie.
2.2 2.2.1
Het productieproces De stappen in het productieproces
Het productieproces bestaat uit vier stappen die hier eerst kort besproken worden. Het proces zoals het hieronder wordt besproken is identiek op alle productielocaties van deze producent.
Breien
Kwaliteitscontrole
Nabewerking
Kwaliteitscontrole
• Parallel • Discontinu • Ongerelateerd • ~ machine • ~ product BOTTLENECK
Relatief beperkte hoeveelheid stock Figuur 2.1: Schematische voorstelling van het breiproces
De eerste stap in het proces is het eigenlijke breien van de stof. Dit gebeurt in breimachines die parallel geschakeld zijn en discontinu produceren. Deze breimachines zijn niet allemaal van hetzelfde type en hebben bijgevolg verschillende productiemogelijkheden. Dit betekent dat niet alle machines alle producten kunnen produceren, en dat deze producten afhankelijk van de machine ook aan een verschillende snelheid wordt geproduceerd. De omschakeling van deze machines tussen twee verschillende producten gebeurt door technici. De tijdsduur van deze omschakeling hangt enkel af van de producten waartussen omgegescha-
Hoofdstuk 2. Probleemsituering
6
keld wordt. Het is belangrijk om hier op te merken dat de tijdsduur van de omschakeling quasi onafhankelijk is van de machines zelf, niettegenstaande verschillende soorten machines worden gebruikt. Deze omschakeltijden kunnen worden benaderd door een onderscheid te maken tussen een omschakkeling binnen eenzelfde productklasse en omschakelingen tussen producten van twee verschillende productklassen. Voor omschakelingen binnen eenzelfde productklasse wordt een omschakelduur berekend op basis van het aantal bobijnen dat moet worden verwisseld. Omschakelingen naar andere productklassen duren heel wat langer dan omschakelingen binnen eenzelfde productklasse, deze omschakeltijd wordt voorlopig met een vaste waarde benaderd. Op de exacte berekening van de omschakeltijden wordt verder ingegaan in sectie 4.2.3. Omdat deze machines een beperkte opslagcapaciteit hebben voor de rollen gebreide stof moeten deze tevens periodisch worden geleegd. Dit kan gebeuren op voorspelbare tijdsintervallen en er wordt verondersteld dat het verhogen van de productietijd met een vast forfait per nodige lediging van de machine volstaat om het proces hier te analyseren. Het is zeer belangrijk om op te merken dat deze ledigingen in de praktijk soms niet onmiddelijk plaatsvinden, machines moeten vaak te lang wachten op beschikbaar personeel voor deze lediging. Alle tijd die machines hierdoor stilstaan is dan ook een verliespost die momenteel zwaar onderschat wordt. De momenten waarop machines geleegd worden zijn deterministisch te voorspellen en kunnen dus door een eenvoudig mechanisme ingepland worden (audiovisuele melding van stilstaande machine). De rest van dit werkstuk veronderstelt dat een dergelijk mechanisme zal worden ge¨ımplementeerd. Omdat dit systeem zo eenvoudig is wordt hier ook niet verder op in gegaan. De rollen die uit de breimachines worden gehaald worden overgebracht naar een kleine tussenvoorrad die zich voor de eerste kwaliteitscontrole bevindt. Vervolgens wordt het product overgebracht naar de tweede stap in het productieproces, de kwaliteitscontrole waar gecontroleerd wordt op fouten in het breiproces. Na de kwaliteitscontrole worden de rollen gestockeerd tot het order volledig geproduceerd is. (Alle stof van eenzelfde order wordt dan verzameld op een dock.) De derde productiestap bestaat uit de chemische nabewerking van de stof. Een zeer grote verscheidenheid aan verschillende nabewerkingen is mogelijk. Voor dit proces worden de orders
Hoofdstuk 2. Probleemsituering
7
meestal gegroepeerd om orders die eenzelfde soort nabewerking vereisen na elkaar te kunnen inplannen met als doel de omschakeltijden te beperken. Deze machine wordt niet uitsluitend gebruikt voor het nabehandelen van producten uit het breiproces en de beschikbaarheid van deze machine is hierdoor moeilijk op voorhand te voorspellen. Historische gegevens tonen dat de wachttijden voor deze machine beperkt zijn, wat betekent dat dit proces voldoende capaciteit heeft om de aangeleverde orders vlot te verwerken. De machine voor de nabewerking van de gebreide stof heeft niet op alle locaties dezelfde productiemogelijkheden. Dit is natuurlijk ook een factor die in overweging moet worden genomen bij het plannen van de orders. In de rest van dit werkstuk wordt er steeds van uitgegaan dat de productiemogelijkheden van de breimachines worden aangepast om ook de productiemogelijkheden van de nabewerkingsmachine op de locatie van de breimachine te weerspiegelen. De vierde en laatste stap van het productieproces bestaat uit de eindcontrole en het verpakken van de orders. De controle gebeurt op dezelfde manier als de controle tussen het breiproces en het nabewerkingsproces, hetzij aan een verhoogde snelheid. Na de controle worden de rollen automatisch naar de verpakkingseenheid getransporteerd en naar het magazijn gebracht tot deze goederen verscheept worden naar de klanten.
2.2.2
Observaties in de productiehal
Geruime tijd werd gespendeerd aan observaties in de productieruimte en interviews met verschillende techniekers om mogelijke oorzaken van vertragingen of ineffici¨enties in het productieproces bloot te leggen. Hieronder worden de meest relevante resultaten kort besproken. De belangrijkste observatie die gemaakt werd is de stilstand van de machines van de eerste productiestap (= het breien zelf). Het betreft hier kapitaalintensieve machines die geruime tijd moeten wachten op een operator die de rol uit de kooi moet halen of die wachten op een omschakeling. Uit gesprekken met operators blijkt dat de productieplanning voornamelijk op ad hoc basis gebeurt. Hieruit kan worden geconcludeerd dat een productieplanning met zowel meer detail (zodat kan worden voorspeld wanneer machines zullen moeten worden geleegd) als een langere tijdspanne (opdat suboptimale keuzes in productopvolgingen, en daarmee onnodig lange omschakeltijden kunnen worden vermeden) voor betere resultaten zou kunnen zorgen.
Hoofdstuk 2. Probleemsituering
8
De hoeveelheid work in progress (WIP) is over het gehele proces zeer beperkt. Tussen de verschillende machines staat telkens een zeer acceptabele hoeveelheid aan producten. Hieruit blijkt duidelijk dat de eerste stap van het productieproces de bottleneck van het proces is (Zie figuur 2.1). De redenering hierachter is de volgende: wanneer alle productiestappen op maximale snelheid produceren zal de productiestap waar de capaciteit het laagste is de beperkende factor vormen voor de output van het volledige proces. Daarom zullen de machines die aan deze productiestap voorafgaan gemiddeld meer output produceren dan deze stap kan verwerken. Dit laatste zal er dan automatisch voor zorgen dat tussenvoorraden zich zullen opstapelen voor deze productiestap. Aangezien zich nergens WIP opstapelt kan er van uitgegaan worden dat de eerste productiestap de bottleneck van het proces is. Zelfs indien just in time (JIT) productiemethodes gebruikt worden, wordt er een voorraadbuffer aangelegd voor de bottleneck die moet vermijden dat de bottleneck niet verder kan produceren wanneer er problemen zijn bij andere machines.
2.2.3
Bepalen van de scope van het probleem
Omdat de eerste productiestap de bottleneck van het proces is kunnen de prestaties van het volledige proces worden verbeterd door de prestaties van de bottleneck te verbeteren (Goldratt & Cox, 2004). De goedkoopste manier om de output van het breiproces te vergroten is het ontwerpen van een nieuwe methode voor het inplannen van orders op deze machines (het alternatief zou het voorzien van meer capaciteit zijn). Door in deze methode rekening te houden met de verschillende eigenschappen van de omschakeltijden, de beschikbaarheid van personeel en het leegmaken van de machines kan de output van het proces worden verhoogd. De voordelen van deze planningsmethode zijn tweeledig. Vooreerst kan de output van het proces worden verhoogd zonder investering in productiecapaciteit en kan aldus beter aan de toenemende vraag worden voldaan. Het tweede voordeel van deze methode is dat de leveringstermijnen veel voorspelbaarder worden dan deze zijn in de huidige situatie. Hierdoor kan worden vermeden dat aan klanten onhaalbare levertermijnen worden beloofd.
Hoofdstuk 2. Probleemsituering
2.3
9
Gedetailleerde beschrijving van het probleem
Figuur 2.2 toont een visualisatie van het probleem dat in dit werkstuk wordt behandeld. Dit schema toont welke jobs op welke machines worden ingepland en hoe veel productietijd nodig is om deze orders te produceren. Daarnaast bevinden de machines zich hier op verschillende geografische locaties, wat zowel een verschil oplevert in afstand naar de klant als een verschil in belasting van het personeel dat de omschakelingen uitvoert. De zwarte blokken duiden de jobs en hun productietijden aan, deze productietijden kunnen verschillen naargelang het type machine waaraan de job wordt toegewezen. De rode zones geven aan hoe veel tijd nodig is om een machine om te schakelen naar het volgende product, deze duurtijden zijn afhankelijk van zowel de voorgaande als de volgende job. time periods
1
2
3
4
5
6
7
8
9
10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30
9001
BE
1
5
9
9002
Machines
BE
10
6
2
8001
RO
3
7
4
8002
11
Staff
RO
Romania Belgium
Available Used Available Used
2 0 2 0
2 0 2 0
2 0 2 0
2 0 2 0
2 0 2 0
8 2 2 2 2
2 2 2 2
2 2 2 0
2 2 2 2
2 0 2 2
2 2 2 2
2 2 2 2
2 0 2 2
2 0 2 2
2 0 2 2
12 2 0 2 2
2 0 2 0
2 4 2 2
2 4 2 2
2 4 2 2
2 2 2 2
2 2 2 0
2 2 2 0
2 0 2 0
2 0 2 0
2 0 2 0
2 0 2 0
2 0 2 0
2 0 2 0
2 0 2 0
Figuur 2.2: Gedetailleerde weergave van schema
De onderste rijen van figuur 2.2 tonen de belasting van de verschillende techniekers die de omschakelingen uitvoeren. Indien deze omschakelingen simultaan ingepland staan zorgt dit ongetwijfeld voor vertragingen van de productie. Zo zien we bijvoorbeeld dat op tijdstippen 18 tot en met 21 twee omschakelingen ingepland staan op machines in Belgi¨e, dit zorgt in dit voorbeeld voor een tekort aan personeel voor deze omschakeling, waardoor de productie van job 9 en/of 2 vertraging zal oplopen. Deze problematiek wordt verder behandeld in hoofdstuk 10 en appendix A.
Hoofdstuk 2. Probleemsituering
10
Het doel van de optimalisatie is hier om de doelfunctie zo laag mogelijk te houden, wat gebeurt als de klanten in verhouding tot hun contributie maximaal tevreden zijn met het tijdig opleveren van hun bestellingen. Dit zorgt zeker indien het probleem op bedrijfsschaal wordt beschouwd voor een zeer complex probleem.
2.4
Huidige planningsmethode
De planning van de breimachines gebeurt momenteel gecentraliseerd maar op een ad hoc basis. Elke dag worden de verschillende openstaande orders beschouwd, wordt er nieuwe informatie verkregen van het commercieel departement en wordt besloten wat die dag geproduceerd zal worden. Voor het schikken van de orders op de machines houdt de planner rekening met een verscheidenheid aan factoren: de productiemogelijkheden van de machines (hierbij worden ook de productiemogelijkheden van de naverwerkingsinstallaties op de verschillende faciliteiten in beschouwing genomen), de beschikbaarheid van grondstoffen, de leveringsdata die beloofd werden aan de klant en het verzekeren van vlotte omschakelingen zijn hierbij de belangrijkste. Het uitbalanceren van deze factoren gebeurt voornamelijk handmatig op basis van ervaring. Dit betekent dat veel van deze informatie ook niet in de vorm van een georganiseerde database beschikbaar is maar op basis van ervaring van de personeelsleden wordt ingeschat. Niettegenstaande deze mensen ervaring hebben met dit proces lijkt het dus onmogelijk dat de uitkomst een globaal optimum kan benaderen. Dit omdat het quasi onmogelijk is om een globaal overzicht te behouden op 300 openstaande orders. Daarnaast geeft de huidige communicatie tussen het commercieel en operationeel departement slechts een uiterst subjectief beeld van het relatieve belang van de orders. Door met behulp van een vaste methode een waarde toe te kennen die het belang van een specifiek order weergeeft zou deze situatie sterk kunnen worden verbeterd, wat zou toelaten om de jobs op een meer objectieve manier in te plannen. Een model dat dit proces optimaliseert kan alle factoren waarmee de planning nu rekening houdt in beschouwing nemen, met daarnaast nog extra factoren.
Hoofdstuk 2. Probleemsituering
2.5
11
Analoge problemen in de literatuur
Deze sectie behandelt de omkadering van het hierboven geformuleerde probleem in de literatuur. Er wordt gestart van een brede kijk op plannen van jobs op machines, waarna meer specifieke bronnen worden behandeld. De bedoeling van deze sectie is niet een exhaustief overzicht te geven van alle literatuur die rond het onderwerp beschikbaar is. Wel wordt het probleem binnen een ruimer kader gesitueerd om vervolgens verschillende oplossingsmethodes van dichtbij te gaan bekijken.
2.5.1
Situering van de aard van het probleem
Een zeer brede omkadering van het probleem dat hierboven wordt omschreven kan gevonden worden in Mirkin & Muchnik (2009). Dit werk beschrijft een multidimensioneel framework waarin planningsproblemen gesitueerd kunnen worden. De verschillende dimensies die in dit werk worden gebruikt om planningsproblemen in te delen worden hieronder kort beschreven en toegepast op het probleem bij deze textielfabrikant. Het eerste onderscheid kan worden gemaakt tussen problemen die precies een productiestap omvatten en problemen die meerdere productiestappen omvatten. (Waarbij problemen met meerdere stappen natuurlijk complexer zijn dan problemen die minder stappen omvatten.) De eigenlijke aard van het breiproces dat in deze paper wordt behandeld omvat meerdere stappen maar voor dit probleem zal enkel de eerste stap - het eigenlijke breiproces - in beschouwing worden genomen. Het probleem wordt hier dus vereenvoudigd naar een singlestage productiesysteem. Een tweede aspect is het aantal machines dat gebruikt wordt in deze productiestap. Hier is het duidelijk dat meerdere machines worden gebruikt, en kan dus gesproken worden over het parallel machine scheduling probleem (PMS). Er kan echter nog een onderscheid gemaakt worden tussen drie verschillende soorten van parallelle machine opstellingen: Identieke parallelle machines Productietijden zijn onafhankelijk van de machine waarop geproduceerd wordt. Uniforme parallelle machines Productietijden zijn afhankelijk van de machine waarop geproduceerd wordt, maar zijn verder identiek (en dus onafhankelijk van het geproduceerde product).
Hoofdstuk 2. Probleemsituering
12
Ongerelateerde parallelle machines Productietijden zijn zowel afhankelijk van de machine als van de job die wordt geproduceerd.
De situatie in het bedrijf is duidelijk van het derde type, de productietijden zijn zowel afhankelijk van de machine waarop geproduceerd wordt als van het product dat wordt geproduceerd. Concreet wordt de maximale productiesnelheid bepaald door het minimum te nemen van het maximale aantal omwentelingen per minuut van de machine en het maximaal aantal omwentelingen per minuut waaraan het product kan/mag worden geproduceerd. Ten derde kan een onderscheid gemaakt worden op het niveau van de jobs. Jobs worden gekarakteriseerd door hun beschikbaarheid (kunnen alle jobs direct worden gestart, of enkel vanaf een bepaald moment in de tijd), hun afhankelijkheid ten opzichte van andere jobs (hiermee wordt bedoeld dat bepaalde jobs afgerond moeten zijn voor een job kan starten) en het al dan niet toestaan van onderbrekingen. In het specifieke geval dat hier wordt behandeld krijgt elke job een specifieke datum voor vrijgave mee, en dient er dus in de planning rekening te worden gehouden dat jobs slechts kunnen starten vanaf dit tijdstip. Hier kan eventueel ook een uitbreiding op het model aansluiting vinden die de beschikbaarheid van de grondstoffen controleert en op basis daarvan de vroegst mogelijke productiedata bepaalt. De jobs zijn echter niet afhankelijk van andere jobs en kunnen dus vrij van de planning van andere jobs worden ingepland. Het onderbreken van de productie wordt uit kwaliteitsoverwegingen eveneens niet toegestaan. (Het opnieuw opstarten van een machine zou kunnen resulteren in een iets andere kwaliteit van de stof.) Ten vierde wordt nog een onderscheid gemaakt tussen online en offline plannen. Bij online plannen wordt de informatie die de planner nodig heeft slechts tijdens het planningsproces zelf beschikbaar, bij offline plannen wordt verondersteld dat alle informatie op voorhand al beschikbaar is. Voor het probleem dat in dit werk wordt behandeld wordt een offline situatie verondersteld, waarbij alle nodige informatie op voorhand beschikbaar is. Doordat de leveringstermijnen vrij lang zijn is dit een aanvaardbare assumptie. Een vijfde en laatste criterium is de gebruikte doelfunctie voor het beoordelen van de jobs. In de literatuur wordt een zeer ruime verscheidenheid aan doelfuncties behandeld (zie sectie 3.1), maar zoals in sectie 3.2 uitgebreid wordt beschreven is hier geopteerd voor een complexere economisch gefundeerde doelfunctie.
Hoofdstuk 2. Probleemsituering
13
Een systematische weergave van de bovenstaande aspecten wordt beschreven door Graham et al. (1979). In dit werk wordt een verkorte notatie beschreven voor de classificatie van planningsproblemen. Naast de bovenstaande indeling is er nog een ander aspect dat nog niet werd behandeld, namelijk de beperkte human resources die beschikbaar zijn voor de omschakeling van de machines. Celano et al. (2008) lossen dit probleem op door gebruik te maken van een simulatie die op basis van eenvoudige beslissingsregels een verwacht resultaat gaat berekenen. (Bijvoorbeeld een first come first serve policy voor de omschakelingen.) Een nadeel van deze aanpak is de gevoelig grotere rekentijden die nodig zijn om een oplossing te evalueren. Dit probleem is verschillend van het probleem zoals het wordt behandeld door Grigoriev et al. (2007), waarbij human resources kunnen worden gebruikt om de productietijden van bepaalde jobs te verkorten. De productiesnelheden in het probleem dat in dit werkstuk wordt behandeld kunnen niet worden versneld door gebruik te maken van extra human resources, maar zijn enkel afhankelijk van de relatie tussen machine en product. De oplossingsstrategie¨en voor problemen van deze aard kunnen in twee grote categorie¨en worden opgedeeld: exacte en (meta-)heuristische technieken. De exacte technieken kunnen met zekerheid een optimale oplossing garanderen, maar zijn vaak te tijdrovend (te lange rekentijd) om op problemen van realistische schaal te worden toegepast. De (meta-)heuristische methoden garanderen geen exacte oplossing maar kunnen wel binnen een redelijk tijdsbestek en ook voor grotere problemen een oplossing van goede kwaliteit leveren.
2.5.2
Exacte oplossingsmethodes
De wiskundige probleemformulering door Tavakkoli-Moghaddam et al. (2009) kan gebruikt worden voor een exacte optimalisatie en ligt aan de basis van het model in hoofdstuk 5 en appendix A. Deze formulering gebruikt een binaire variabele van de vorm xjkm , die als betekenis heeft dat wanneer deze binaire variabele de waarde 1 heeft job j ingepland is op positie k van machine m. Een alternatieve exacte modelformulering kan gevonden worden in Gharehgozli et al. (2009). Hier wordt de voornaamste beslissingsvariabele niet gespecifieerd in functie van de positie, maar afhankelijk van de voorgaande job. De binaire variabele xm ij is hier gelijk aan 1 wanneer job j net na job i gepland staat op machine m.
Hoofdstuk 2. Probleemsituering
2.5.3
14
(Meta-)Heuristische oplossingsmethodes
Uit de literatuur die (meta-)heuristische oplossingsmethodes beschrijft is duidelijk dat er een enorm grote verscheidenheid aan oplossingsmethodes bestaat voor dit probleem. Hier wordt niet geprobeerd een exhaustief overzicht van al deze technieken weer te geven, maar worden enkele recente papers rond dit onderwerp aangehaald. Een duidelijk onderscheid kan worden gemaakt tussen twee grote categorie¨en van metaheuristieken: heuristieken gebaseerd op lokaal zoeken en heuristieken gebaseerd op populaties. Bij lokaal zoeken start het algoritme van een bepaalde oplossing (die op verscheidene manier kan worden samengesteld) waarin relatief kleine wijzigingen worden aangebracht. Vervolgens wordt de licht gewijzigde oplossing ge¨evalueerd om te zien of deze wijziging al dan niet voor een betere evaluatie van het totale schema gaat zorgen. De meest voorkomende technieken die gebruik maken van lokaal zoeken zijn: Simulated annealing Threshold accepting Tabu Search Variable neighborhood search
Daartegenover staat het populatie gebaseerd zoeken, waarbij wordt gestart van een groot aantal verschillende schemas (die tevens op zeer uiteenlopende manieren kunnen worden gegenereerd). Door het combineren van deze verschillende oplossingen - waarbij voorrang wordt gegeven aan oplossingen met betere evaluaties - worden nieuwe schema’s gemaakt. De meest voorkomende technieken die gebruik maken van populaties zijn: Scatter search Genetic algorithm Ant colony optimization
Enkele recente voorbeelden van auteurs die metaheuristieken gebaseerd op lokaal zoeken gebruiken om het PMS probleem op te lossen zijn: Bilge et al. (2004) (Tabu search), Chen
Hoofdstuk 2. Probleemsituering
15
& Wu (2006) (Threshold accepting), Yazdani et al. (2010) (Variable neighborhood search), Celano et al. (2008) en Kim et al. (2002) (Simulated Annealing). De paper van Anghinolfi & Paolucci (2007) behandelt een hybride techniek gebaseerd op verschillende lokaal zoekende metaheuristieken. Naast de uitgebreide verzameling aan werken die de PMS problematiek aanpakt met op lokaal zoeken gebaseerde technieken is er ook een zeer uitgebreide verzameling aan auteurs die dit probleem met behulp van populatie-gebaseerde metaheuristieken oplossen. Een voorbeeld hiervan is Kumar et al. (2011) waarin een zeer eenvoudige ant colony optimization wordt beschreven voor het oplossen van het PMS probleem. De meeste literatuur focust zich echter op genetische algoritmen en varianten hierop. Een recent werk is dat van Vallada & Ruiz (2011) waarin de ontwikkeling en uitbreidingen op een genetisch algoritme voor het PMS probleem uitgebreid worden besproken. Een van de interessante uitbreidingen is gebaseerd op Kurz & Askin (2001) waarin lokaal zoeken wordt gebruikt om de kwaliteit van de nakomelingen te verbeteren. Daarnaast wordt naar analogie met Ruiz & Allahverdi (2007) de n-tournament methode gebruikt om de ouders te selecteren, een techniek die ook volgens Glover & Kochenberger (2003) betere resultaten geeft dan de klassieke roulette selectie mechanismen. Een hybride versie tussen beide categorie¨en wordt ge¨ıllustreerd in Carvalho & Ceccarelli (2001). Hier wordt een memetic algoritme gebruikt dat zowel populatie als lokaal zoeken gebruikt om een oplossing te zoeken.
2.5.4
Conclusies uit literatuuronderzoek
Afgaand op het bovenstaand literatuuronderzoek kan worden besloten dat zeer uiteenlopende technieken mogelijk zijn om dit probleem te gaan oplossen. Om een duidelijk beeld te krijgen van de haalbaarheid en de performantie van de verschillende technieken komen in deel II enkele verschillende oplossingsstrategie¨en aan bod. De eerste oplossingsmethode is de exacte oplossingsmethode (zie hoofdstuk 5). Deze wordt bekomen door middel van lineair programmeren met als basis de werken die hierboven vermeld werden. Deze modellen worden vervolgens opgelost met behulp van Gurobi 4.5.1. Met betrekking tot de heuristische oplossingsmethodes wordt gekozen om zowel een heuristiek gebaseerd op lokaal zoeken (hoofdstuk 8) als een heuristiek op basis van populaties (hoofd-
Hoofdstuk 2. Probleemsituering
16
stuk 7) te implementeren. Celano et al. (2008) argumenteert dat problemen waarbij telkens relatief zware berekeningen nodig zijn voor het berekenen van de doelfunctie soms beter met metaheuristieken gebaseerd op lokaal zoeken worden uitgerekend omdat populatie gebaseerde technieken dan een zeer grote rekentijd hebben. Voor beide gevallen wordt gekozen voor traditionele technieken. Voor populatie gebaseerde technieken is dit het genetisch algoritme, voor de lokaal zoeken categorie is dit simulated annealing.
Hoofdstuk 3
Doelfunctie en inputparameters Na de identificatie van het probleem in het vorige hoofstuk wordt nu besproken hoe de oplossingen van het probleem ge¨evalueerd zullen worden en hoe dit bedrijfsspecifieke probleem kan worden omgezet naar een formulering met meer generieke parameters. Dit hoofdstuk bestaat uit drie delen. Het eerste deel omschrijft de meest traditionele doelfuncties uit de literatuur. In het tweede deel wordt een logische formulering gegeven van de doelfunctie die hier zal worden gebruikt. Vervolgens wordt in het derde en laatste deel deze logische formulering in een wiskundige formulering omgezet. Deze doelfunctie werd opgesteld met de medewerking van het commercieel departement van het textielbedrijf.
3.1
Frequent gebruikte doelfuncties uit de literatuur
Voor de klassieke formulering van het PMS probleem worden verschillende soorten doelfuncties gebruikt, die meestal een minimalisatie van een (groep van) variabele(n) inhouden. Een overzicht hiervan wordt gegeven door Pfund et al. (2008) waarin twee grote categorie¨en van doelfuncties worden onderscheden: enkelvoudige objectieven en meervoudige objectieven. Hieronder een kort overzicht van de enkelvoudige objectieven die in dit werk worden aangehaald. (Deze lijst is niet exhaustief, maar geeft een goed overzicht van de meest gebruikte objectieven.) De onderstaande objectieven worden allemaal geminimaliseerd. Makespan De totale tijdsduur die nodig is voor het produceren van alle jobs. (= het tijdstip waarop de laatste job afgewerkt is.) 17
Hoofdstuk 3. Doelfunctie en inputparameters
18
Sum of the weighted completion times Dit objectief laat toe om verschillende jobs een andere weging te geven, waardoor belangrijkere jobs zwaarder doorwegen in de oplossing. Maximum tardiness De parameter “tardiness” wordt hier gespecifieerd als het maximum van nul en het verschil van de beloofde datum en de datum van oplevering voor alle jobs. Weighted tardiness Hier wordt ervoor gezorgd dat de tardiness gebalanceerd wordt over alle jobs, waarbij belangrijkere jobs opnieuw zwaarder gaan doorwegen. Number of tardy jobs (= Lateness) Het aantal jobs dat later dan de beloofde datum wordt opgeleverd. Weighted number of tardy jobs Analoog aan de voorgaande, maar met grotere gewichten voor belangrijkere jobs.
Situaties in de praktijk hebben echter vaak meervoudige objectieven nodig omdat niet alle doelstellingen in een enkel objectief kunnen worden samengevat. Om de werkelijke doelstellingen van het bedrijf te gaan weerspiegelen wordt dan vaak een gecombineerde doelfunctie ontwikkeld die enkele (zo min mogelijk conflicterende) objectieven gaat combineren. Een goed opgestelde doelfunctie moet ervoor gaan zorgen dat de output van het model quasi exact kan worden aangenomen zonder manuele aanpassingen. Met deze redenering in het achterhoofd wordt voor dit probleem een formulering ontworpen die alle voor het bedrijf relevante informatie probeert op te nemen. De doelfunctie die wordt ontwikkeld is in essentie een combinatie van de “weighted tardiness” en de “weighted lateness” enkelvoudige objectieven die hierboven omschreven werden.
3.2
Logische Samenstelling van de doelfunctie
Het optimale productieschema wordt gezocht door het minimaliseren van een kostenfunctie die de kosten gemaakt door het laattijdig leveren van orders weergeeft. Deze doelfunctie bestaat uit een vaste en een variabele component. Hieronder wordt verstaan een kost die wordt gemaakt vanaf het ogenblik dat een order te laat wordt opgeleverd, respectievelijk de kost die wordt gemaakt per dag dat het order te laat wordt opgeleverd.
Hoofdstuk 3. Doelfunctie en inputparameters
19
Totale kosten laattijdige orders
Variabele kosten
Commercieel risico
Kans op contributieverlies
Benodigde werkuren
Administratief
Loonkost VTE
Omvang van order
Afstand naar klant
Transport
ABC score (contributie)
Vaste kosten
Order-afhankelijk Order-onafhankelijk
Figuur 3.1: Samenstelling van de vaste en variabele kosten in de doelfunctie
De vaste component kan verder worden opgedeeld in de transportkost en de administratieve kosten die worden gemaakt. Wanneer een order te laat wordt afgewerkt moet het bedrijf zelf instaan voor de spoedlevering naar de klant van de beloofde producten. Administratieve kosten bestaan uit zowel interne als externe communicatie bij laattijdig leveren van orders door personeelsleden die nodig is wanneer een order te laat is. Transportkosten zijn afhankelijk van zowel de afstand naar de klant als de omvang van de bestelling. Concreet komt dit tot uiting in een bepaalde kost per kilometer die afhankelijk is van het gewicht van de order. Voor orders tot 400 kilogram is dit 0,38 EUR/km, orders tot 800 kilogram is dit 0,45 EUR/km en voor orders tot 1500 kilogram bedraagt dit 0,48 EUR/km. Voor de berekening van de administratieve kosten wordt uitgegaan van de gemiddelde kost voor een voltijds equivalent van 55.392 EUR op jaarbasis (deze informatie werd verkregen van accounting). Wanneer dit wordt herrekend naar loonkost per nuttig werkuur wordt een bedrag van 31,65 EUR bekomen (uitgaand van 250 gewerkte dagen per jaar en 7 nuttige werkuren per gewerkte dag). Gebaseerd op interviews met personeelsleden kan de gemiddelde duur van externe communicatie op dertig minuten worden vastgelegd. Voor interne communicatie bedraagt dit eveneens dertig minuten, maar daar dit steeds twee personeelsleden betreft wordt deze tijdsduur dubbel geteld voor het bepalen van de kosten. De totale vaste administratieve kost wordt aldus vastgelegd op 47,48 EUR (= 1,5 nuttige werkuren).
Hoofdstuk 3. Doelfunctie en inputparameters
20
De variabele component bestaat uit het commercieel risico dat het bedrijf loopt door te laat te leveren bij klanten. Op basis van een ABC-categorisatie van de klanten werd de gemiddelde contributie van deze klanten bepaald. Op basis van interviews met het commercieel departement van het bedrijf werd bepaald dat een levering een maand na beloofde datum in twintig procent van de gevallen zal leiden tot een verlies van de referentie, een verlies van volume komt eveneens in twintig procent van de gevallen voor en in twintig procent van de gevallen zal de verkoper in de toekomst zijn marge moeten verlagen. Zowel margeverlies als volumeverlies komt meestal neer op 10% contributieverlies voor de specifieke klant. De tijdshorizon voor deze berekening wordt vastgelegd op een jaar om geen overschatting te maken van de kosten. Gemiddelde contributies voor A, B en C klanten bedragen op jaarbasis respectievelijk 90.000, 10.000 en 225 EUR (deze getallen werden berekend op basis van gedetailleerde omzetcijfers van 2010). Wanneer twintig werkdagen per maand worden verondersteld en we de kosten op maandbasis over deze verdelen bekomen we een kost van 1.000, 115 en 3 EUR per dag laattijdige levering voor A, B en C klanten respectievelijk. Voor implementatie in de praktijk kan dit principe verder worden uitgewerkt met meer gedetailleerde wegingen. Zo kunnen de echte contributies van de klanten als parameter worden gebruikt en niet de gemiddelden per categorie. Daarnaast kunnen ook grotere waarden worden voorzien voor klanten uit risicocategorie¨en zoals dochterondernemingen van klanten met een grote contributie of nieuwe klanten die in de toekomst potentieel grote orders zullen plaatsen.
3.3
Wiskundige formulering van de doelfunctie
De bovenstaande redenering kan worden weergegeven in de onderstaande niet-lineaire vergelijking 3.1. Deze formule is aangepast voor gebruik in het exacte optimalisatiemodel van hoofdstuk 5. (Noot: in hoofdstuk 5 wordt deze formule herschreven in lineaire vorm, maar voor de eenvoud wordt hier eerst de eenvoudiger te interpreteren niet-lineaire vorm gepresenteerd.) NJ N M X NK X X [Fjm ∗ xjkm ] ∗ lj ) + (Vj ∗ tj )] M in. [( j=0 m=0 k=0
Waarvoor geldt dat:
(3.1)
Hoofdstuk 3. Doelfunctie en inputparameters
21
Fjm De vaste kost gemaakt wanneer de job j die op machine m ingepland staat te laat wordt opgeleverd. xjkm Beslissingsvariabele die gelijk is aan 1 wanneer job j op de k-de positie ingepland staat op machine m. lj Binaire variabele die gelijk is aan 1 wanneer job j te laat wordt opgeleverd. Vj Variabele kost gassocieerd met job j. tj Integer variabele die aangeeft hoe veel dagen job j te laat werd opgeleverd.
Het eerste deel van deze vergelijking (
N M N K P P
[Fjm ∗ xjkm ] ∗ lj ) geeft de vaste kosten weer die
m=0 k=0
gemaakt worden. Hierbij wordt eerst bepaald welke vaste kost van toepassing is, omdat deze afhankelijk is van de machine waarop de job staat ingepland. De reden hiervoor is dat de productielocatie wordt bepaald aan de hand van de machine waarop de job staat ingepland, deze locatie heeft dan op zijn beurt een invloed op de transportkost (zie figuur 3.1). Vervolgens wordt deze vaste kost vermenigvuldigd met de variabele die aangeeft of de job al dan niet te laat werd opgeleverd. Het tweede deel van de vergelijking (Vj ∗ tj ) geeft de variabele kosten weer. Deze kosten zij afhankelijk van de job en worden vermenigvuldigd met een variabele die aangeeft hoe veel dagen een bepaalde job te laat wordt opgeleverd. Belangrijk voor de interpretatie van deze en alle volgende modelformuleringen is dat alle tijdstippen waarover wordt gesproken in dagen worden uitgedrukt. Dit gaat zowel over productietijden, beloofde levertijden, release times en alle andere tijdstippen. Enkel in hoofdstuk 10 en bij het model in appendix A en bij het detecteren van conflicten in hoofdstuk 10 worden tijdsvariabelen gebruikt die periodes van een kwartier aanduiden, deze worden gebruikt voor het opsporen van personeelsconflicten bij het omschakelen van machines.
Hoofdstuk 4
Data generator De kwaliteit van de oplossing hangt sterk af van de input die gebruikt wordt om de oplossingsmethode te testen. Daarom werd hier geopteerd om een volledig onafhankelijke data generator te ontwikkelen die aan de hand van verschillende parameters willekeurige, maar representatieve datasets kan aanleveren. De data nodig voor de verschillende modellen wordt bondig samengevat in tabel 5.5. (En voor de uitgebreide versie in tabel A.5.) Vaste kosten Variabele kosten
Tekort productietijd Productietijd
Historische informatie
Personeelscode
Productiemogelijkheden
Release time
Productiesnelheden (Rev/min)
Due date
Batchgrootte Producttypes Productgroepen
Data Generator
Systeemeigenschappen
Beschikbaarheid personeel Machine locatie Nationaliteiten
ABC Scores Nationaliteiten
ABC Score
Transportafstanden
Ordergrootte
Ordergroottes
Parameter Object
User
# Jobs # Machines
Machine type
Machinegroepen # Cm / omwenteling
Product type
Figuur 4.1: Schematische weergave van input en output van de data generator
In dit hoofdstuk wordt de werking van deze data generator kort besproken, zonder in technisch detail te gaan over de implementatie van deze generator in C++. Dit gebeurt door eerst kort de bedrijfsdata te bespreken waarop de generator gebaseerd is, om dan vervolgens de sequenti¨ele stappen in de generatorprocedure kort te omschrijven. 22
Hoofdstuk 4. Data generator
4.1
23
Input voor data generator
4.1.1
Namen voor illustratieve doeleinden
Met als doel de uitkomst van de verschillende optimalisaties met een groter realisme te kunnen weergeven worden de namen van machines en productkwaliteiten ge¨ımporteerd. Deze worden gekoppeld aan gegenereerde nummers die dan gedurende de rest van de uitvoering gebruikt worden, wanneer output gegeneerd is kunnen de namen opnieuw aan de output kunnen worden gekoppeld.
4.1.2
Productiespecificaties
De relaties tussen de verschillende productkwaliteiten en de machines wordt vervat in drie dimensies. Productiemogelijkheden: Welke producten kunnen door welke machines al dan niet worden geproduceerd? Maximale productiesnelheid: Wat is het grootste aantal omwentelingen per minuut dat kan worden gehaald voor een bepaald product op een bepaalde machine. Maximale en minimale rolgrootte: De minimale rolgrootte is het kleinste aantal omwentelingen dat moet gedraaid worden alvorens de machine kan worden stilgelegd. Dit is vooral van belang voor het bepalen van de eerste kwaliteitscontrole die gebeurt na het produceren van een minimaal aantal meter. De maximale rolgrootte is van belang voor het bepalen van het aantal keer dat de machine leeggemaakt zal moeten worden.
4.1.3
Categorie¨ en van kwaliteiten
De omschakelingen die gebeuren kunnen klein of groot zijn. Een kleine omschakeling is vooral afhankelijk van het aantal bobijnen dat moet worden gewisseld, maar een grote omschakeling is hiervan niet langer afhankelijk en wordt een vaste productietijd toegewezen. (Deze bedraagt ongeveer een volledige shift van acht uur). Om te bepalen wanneer een grote omschakeling plaatsvindt en wanneer een kleine omschakeling plaatsvindt worden de kwaliteiten ingedeeld in veertien verschillende productcategorie¨en. Wanneer een omschakeling plaatsvindt tussen producten uit verschillende categorie¨en, wordt deze beschouwd als een grote omschakeling.
Hoofdstuk 4. Data generator
24
Deze categorie¨en zijn gecree¨erd in functie van dit onderzoek, voordien bestond geen enkele data die de omschakeltijden op een systematische manier kon gaan voorspellen. Alle planning gebeurt vandaag dus op gevoel en ervaring, wat natuurlijk geen optimale oplossingen kan garanderen en daarinboven een zware werklast vormt voor de werknemers. Voor de implementatie is het belangrijk om op te merken dat het bovenstaande een vereenvoudiging is van de echte situatie. Voor een implementatie van het systeem zou het (op termijn) noodzakelijk zijn om een meer uitgebreide studie te maken van de factoren die de omsteltijden bepalen.
4.1.4
Instellingen van systemen
Een breimachine beschikt over veertien systemen waarop draden worden aangesloten. Twee van deze systemen zijn specifiek bedoeld voor vuldraden. De instellingen die voor de berekeningen van belang zijn, zijn de verschillende soorten bobijnen en vuldraden die op de verschillende systemen worden geplaatst om een specifieke kwaliteit te produceren. Hiervan wordt zowel het aantal als het type ge¨ımporteerd. Deze informatie is nodig om bij kleine omschakelingen (binnen eenzelfde categorie van producten) de duur van de omschakeling te gaan bepalen.
4.1.5
Nationaliteiten en ABC scores
Om voor de verschillende orders een nationaliteit en een ABC score te kunnen genereren werd ook informatie verzameld over de verdeling van de omzet en contributiecijfers over verschillende landen en ABC scores heen. De exacte interpretatie en verwerking van deze gegevens komt later aan bod.
4.1.6
Transportafstanden
Voor de verschillende productie en klantlocaties werden proxies berekend die de afstand naar de klant gaan voorstellen. Voor deze oefening is het onnodig om deze informatie al te precies te maken, maar voor implementatie binnen het bedrijf zal het natuurlijk aan te raden zijn om deze informatie te verfijnen.
Hoofdstuk 4. Data generator
4.1.7
25
Ordergroottes
Informatie over de omvang van de orders (in aantal meter) voor de verschillende klanten werd eveneens ge¨ımporteerd in de data generator. Op basis van deze informatie kan dan een realistische verdeling van de ordergroottes worden bepaald. Verder wordt binnen de productieafdeling ook vooropgesteld dat de maximale ordergrootte om en bij de 2000 meter bedraagt. Deze opsplitsing gebeurt al bij de interne verwerking van de orders. De reden hiervoor is dat grotere orders op verschillende machines geproduceerd moeten worden (om zeer lange productietijden te vermijden), wat kleine verschillen in kwaliteit tot gevolg kan hebben.
4.1.8
Machinegroepen
Het aantal aanwezige breimachines van verschillende types wordt ook ge¨ımporteerd, met als doel een realistische verdeling van verschillende machines te hebben voor grotere productieproblemen.
4.1.9
Centimeters per omwenteling
Een andere belangrijke eigenschap van de verschillende producten is het aantal centimeter stof dat per omwenteling wordt geproduceerd. Deze informatie wordt samen met de informatie over de machine snelheid en batchgrootte gebruikt voor het bepalen van de productieduurtijden.
4.2
Procedure voor datageneratie
Om willekeurige en representatieve informatie te kunnen genereren werd een functie ontwikkeld die op basis van enkele inputs een random dataset gaat genereren. De inputs die hierbij gebruikt worden zijn het aantal jobs, het aantal machines, de beschikbare tijd om de orders te produceren en een staffing code. De beschikbare tijd is een parameter die wordt gebruikt om te bepalen hoe veel tijd er tekort is om alle jobs op het beloofde tijdstip af te werken. Dit wordt berekend aan de hand van de minimale tijd die nodig is om alle orders te produceren in een fictieve optimale situatie. Deze situatie wordt bekomen door voor alle jobs de minimale productietijden te nemen, te verhogen met een minimale omschakeling en deze totale som vervolgens te delen door het
Hoofdstuk 4. Data generator
26
aantal aanwezige machines. Deze minimale tijd kan dan vervolgens naar wens worden vergoot of verkleind om situatie met meer of minder druk op het productieapparaat weer te geven. Er kan ook een staffing code worden meegegeven, aan de hand van deze code wordt een bepaald patroon van beschikbaarheid voor omschakelend productiepersoneel gekozen. dit kan bijvoorbeeld een continue aanwezigheid van 1 persoon zijn, of de aanwezigheden van 2 personen overdag en 1 persoon ’s nachts.
4.2.1
Informatie over orders
De eerste stap van de procedure is het toekennen van nationaliteiten aan de verschillende orders die in de dataset zullen worden opgenomen. De toekenning gebeurt op basis van strata vastgelegd op basis van omzetcijfers uit 2010. Op deze manier wordt ook voor kleinere datasets een verdeling bekomen die sterk overeenkomt met de verdeling in realiteit. Vervolgens worden voor de verschillende orders ABC scores toegekend. Deze ABC scores duiden het commercieel belang van de verschillende klanten aan. Deze scores worden toegekend op basis van willekeurig gegenereerde nummers die worden vergeleken met het relatief aantal klanten in de verschillende categorie¨en voor de verschillende landen. (Zo heeft een Franse order relatief meer kans om aan een A klant toe te behoren dan een order uit Roemeni¨e.) Nadat alle ABC scores bepaald zijn worden de orders een ordergrootte toegewezen. Dit gebeurt op basis van de ABC scores. Hier komt het er op neer dat A klanten grotere orders plaatsen dan B of C klanten. Deze toekenning gebeurt opnieuw op basis van strata om ook voor kleine datasets een zo representatief mogelijke verdeling na te streven. Voor alle orders wordt eveneens een product gekozen. Deze keuze gebeurt willekeurig uit de groep van producten die geproduceerd kan worden door de machines opgenomen in deze dataset en waarvan eveneens voldoende informatie over de bobijninstellingen beschikbaar zijn. Dit laatste is een belangrijk element voor de implementatie van dit project, het is van groot belang dat de informatie ook voor alle varianten op producten correct wordt aangevuld. Voor het bepalen van de beloofde tijdstippen van opleveren en de tijdstippen waarop orders in het systeem worden ge¨ıntroduceerd wordt er gebruik gemaakt van de “minimal makespan” waarover hierboven al meer uitleg gegeven werd.
Hoofdstuk 4. Data generator
27
Om tijdstippen vast te leggen voor de release dates en de due dates van de orders worden dan willekeurig tijdstippen gekozen in de eerste respectievelijk tweede helft van de gereduceerde makespan. Eventuele uitbreidingen van het model die rekening houden met de beschikbare garens kunnen hier hun aansluiting vinden door de release times aan te passen naargelang de beschikbaarheid van de garens.
4.2.2
Machine informatie
Voor de machines waarmee geproduceerd wordt zijn er twee aspecten die van belang zijn, vooreerst het type van de machine, wat invloed heeft op de productiemogelijkheden en productiesnelheden voor verschillende productkwaliteiten, en de locatie van de machine, die invloed heeft op de transportkost van de expresverzendingen voor laattijdig geleverde orders. Het type van de machine wordt willekeurig bepaald, met kansen die gebaseerd zijn op het aantal machines dat in realiteit aanwezig op alle verschillende productielocaties. De productielocaties van de machines worden op hun beurt bepaald door een vooraf vastgelegde verhouding die aanduidt hoeveel van de machines zich in Roemeni¨e en hoeveel machines zich in Belgi¨e bevinden. Voor de eenvoud werden in deze studie slechts twee productielocaties opgenomen. Het is echter perfect mogelijk om meer dan twee productielocaties te gebruiken zonder structurele aanpassingen aan de oplossingsmodellen.
4.2.3
Parameters met betrekking tot productie
De twee types parameters die betrekking hebben tot de productie zelf zijn de productietijden en de omschakeltijden. Dit zijn tevens de twee parameters die het meest complex zijn om te bepalen. De productietijd van de jobs bestaat uit twee componenten: de tijd die de machine nodig heeft om te produceren en de tijd die nodig is om de rollen uit de machine te verwijderen wanneer deze vol zijn. Als eerste stap voor het bepalen van de productietijden worden de productiesnelheden omgezet van aantal omwentelingen per minuut naar aantal meters per minuut. Hiervoor wordt
Hoofdstuk 4. Data generator
28
het minimum genomen van de maximale snelheid van de machine in aantal omwentelingen per minuut en de maximale productiesnelheid voor het specifieke producttype (eveneens in omwentelingen per minuut). Dit getal wordt dan vermenigvuldigd met het aantal meter stof dat de order groot is om de totale netto productietijd te bekomen. De tweede stap is dan het bepalen van het aantal keren dat de machine wordt stilgelegd om de rollen uit de machine te verwijderen en het bepalen van de tijd die nodig is om de rollen uit de machine te gaan verwijderen. Machines worden een eerste maal gestopt na de minimale batchgrootte om de kwaliteit te controleren, deze controle neemt een tiental minuten in beslag. Vervolgens wordt de machine telkens gestopt wanneer de maximale batchgrootte bereikt is. Voor deze handelingen wordt telkens rekening gehouden met een vijftal minuten extra productietijd alvorens de machine opnieuw kan produceren. Belangrijk is hier dat in het huidige productieregime het ledigen van de machines niet als prioritair wordt aanzien. Om maximale productiviteit te kunnen bekomen zou dit echter wel zo moeten zijn. Verschillende observaties in de productiehal hebben aangetoond dat machines soms lang moeten wachten om geledigd te worden, ondanks dat deze handeling weinig tijd in beslag neemt en daarenboven perfect voorspelbaar is. Om te bekomen dat het ledigen van deze machines correct en tijdig gebeurt zijn nieuwe controlemechanismen nodig waarbij specifieke personen voor specifieke machines worden verantwoordelijk gesteld. Daarbij kan nog een systeem komen die een signaal geeft wanneer een machine gestopt is met produceren omdat de kooi geleegd dient te worden. Zeer grote en zeer eenvoudig te bekomen productiviteitswinsten kunnen op deze manier worden bekomen. Deze twee componenten worden vervolgens samengevoegd en omgezet naar productietijden in dagen (dit omdat alle andere tijdstippen in het model in aantal dagen worden weergegeven). De omschakeltijden tussen de verschillende jobs worden twee aan twee berekend. De eerste stap in deze berekening is controleren of de twee kwaliteiten al dan niet tot eenzelfde categorie behoren. Indien deze tot verschillende categorie¨en behoren is er automatisch sprake van een grote omschakeling. (Een grote omschakeling is een omschakeling van acht uur.) Indien de twee kwaliteiten tot eenzelfde categorie behoren wordt de belangrijkste driver voor de duur van de omschakeling bepaald, dit is het aantal bobijnen dat verwisseld moet worden. Deze informatie wordt gehaald uit de informatie van de tien verschillende systemen en
Hoofdstuk 4. Data generator
29
vuldraden aanwezig op de machine. Elk van deze tien systemen kan een aantal bobijnen van een bepaald type hebben. Per systeem wordt nagegaan of de bobijnen van eenzelfde of een verschillend type zijn. Indien de bobijnen van een verschillend type zijn wordt de grootste van de twee aantallen genomen als het te verwisselen aantal bobijnen. Wanneer de bobijnen echter van eenzelfde type zijn wordt er enkel rekening gehouden met het verschil in aantal bobijnen op het systeem. Om dan voor het geval waarbij de kwaliteiten tot eenzelfde categorie behoren eveneens tot een bepaald aantal tijdsperiodes te komen wordt het aantal bobijnen vermenigvuldigd met het aantal minuten nodig voor het wisselen van een bobijn en gesommeerd voor alle bobijnen (0,8 minuten / bobijn). Dit wordt dan verhoogd met een constante forfait die de andere handelingen voor een omschakeling voorstelt. Vervolgens worden deze waarden omgezet naar aantal tijdsperioden van een kwartier (deze worden naar boven afgerond). Deze tijdsperioden kunnen waar nodig worden omgerekend naar duurtijden in aantal dagen. Deze berekening wordt ge¨ıllustreerd in tabel 4.1. Deze tabel toont de berekening voor de omschakeling van producttype BAMBOO-2C 24 naar producttype APOLLO-C 14. De rijen in deze tabel geven de tien verschillende systemen weer. In de tweede en derde kolom zien we dan de invulling van deze systemen voor het producttype BAMBOO-2C 24. Daarnaast wordt opnieuw dezelfde informatie herhaald voor APOLLO-C 14. Op het eerste systeem zijn de bobijnen van een verschillend type, wat betekent dat alle bobijnen van dit systeem zullen moeten worden verwijderd en dertig nieuwe bobijnen geplaatst zullen moeten worden. Let wel dat hier enkel de bobijnen die worden opgezet worden verrekend voor het schatten van de omschakeltijd.
4.2.4
Vaste kosten voor laattijdig leveren
Het bepalen van de vaste kost voor laattijdige levering van de verschillende orders is de volgende stap in het generatieproces. Zoals in ge¨ıllustreerd wordt in figuur 3.1 bestaan de vaste kosten uit twee componenten: de transportkost en de administratieve kost. Voor het bepalen van de transportkost wordt eerst bepaald of het om een klein, middelgrote of grote order gaat. Afhankelijk van de categorie wordt er dan een kost per kilometer geselecteerd die dan voor alle machines wordt vermenigvuldigd met de afstand die moet worden overbrugd voor de levering. Deze afstanden verschillen van machine tot machine omdat deze zich op verschillende productielocaties bevinden.
Hoofdstuk 4. Data generator System
30
BAMBOO-2C 24
APOLLO-C 14
Same Y/N
To change
Type
# Bobins
Type
# Bobins
1
167504-1
15
193004-90
30
No
30
2
0
0
0
0
Yes
0
3
182019-1
15
790759-1
15
No
15
4
60059-1
15
60059-1
15
Yes
0
5
182019-1
15
790759-1
15
No
15
6
791504-32
15
791504-32
15
Yes
0
7
0
0
0
0
Yes
0
8
0
0
0
0
Yes
0
9
0
0
0
0
Yes
0
10
0
0
0
0
Yes
0
TOTAL
60
Tabel 4.1: Berekening van omschakeltijden tussen twee kwaliteiten binnen eenzelfde categorie.
Nadat de transportcomponent van de vaste kosten bepaald is wordt deze nog verhoogd met de administratieve component die in sectie 3.2 werd vastgelegd op 47,48 EUR. (= de bedrijfskost voor 1,5 nuttige werkuren). Op deze manier wordt voor elk order en elke machine een bepaalde vaste kost bepaald voor het laattijdig leveren.
4.2.5
Variabele kosten voor laattijdig leveren
De variabele kosten voor laattijdige levering worden op eenvoudige manier bepaald op basis van de ABC scores. Voor de verschillende ABC scores werd een bepaalde commerci¨ele kost berekend per dag laattijdig leveren. (Deze kosten worden eveneens berekend in sectie 3.2) Alle orders worden eenvoudigweg toegewezen aan een van deze drie categorie¨en voor variabele kosten.
4.2.6
Maximale waarden voor indices
De maximale waarden voor de indices die de jobs, machines en productielocaties voorstellen (NJ, NM en NP) worden automatisch geven wanneer een groep parameters wordt aangemaakt. Voor het aantal posities op de verschillende machines (NK) wordt de waarde van NJ gebruikt,
Hoofdstuk 4. Data generator
31
dit zou de situatie zijn waarbij alle jobs op eenzelfde machine worden ingepland, en dit is dan ook meteen het worst case scenario voor deze index. Voor het aantal tijdsperiodes (NT) is de berekening echter ietwat complexer. (Deze parameter wordt enkel gebruikt voor het uitgebreide model in appendix A) De eerste stap in de berekening is het selecteren van alle maximale productietijden en deze sommeren. Vervolgens wordt er een buffer toegevoegd voor alle omschakeltijden, voor de eenvoud wordt hier verondersteld dat het hier telkens om grote omschakelingen gaat. Deze tijden worden vervolgens toegevoegd aan de som van de omschakeltijden. Een derde en laatste component van deze som is een buffer voor mogelijke wachttijden op beschikbaar personeel. Hier wordt opnieuw gekozen voor een worst case scenario waarbij alle omschakelingen gelijktijdig dienen te gebeuren door een enkel personeelslid. (Dit betekent dat het eerste order geen wachttijd heeft, het tweede order een wachttijd meekrijgt gelijk aan een grote omschakeling, het derde order een wachttijd meekrijgt gelijk aan twee grote omschakelingen, etc.) Deze totale wachttijd wordt aan de som toegevoegd en deze wordt herrekend naar tijdsperiodes van een kwartier, deze waarde wordt dan gebruikt als waarde voor NT.
Deel II
Ontwerpen van oplossingsmethodes
32
Hoofdstuk 5
Exact optimalisatiemodel Dit hoofdstuk bespreekt de gebruikte exacte optimalisatiemethode. De modelformulering is gebaseerd op het werk van Tavakkoli-Moghaddam et al. (2009). Dit hoofdstuk start met de wiskundige formulering en presenteert daarna kort de resultaten voor verschillende probleemgroottes. Dit model werd opgelost met Gurobi 4.5.1 gebruikmakend van de C++ interface. De duiding bij het model wordt in dit hoofdstuk tot een minimum beperkt om de leesbaarheid te verbeteren. Voor een uitgebreide uitleg wordt doorverwezen naar hoofdstuk 6. Naast het model dat in dit en het volgende hoofdstuk wordt beschreven werd een complexer model ontwikkeld dat tevens rekening houdt met mogelijke personeelstekorten voor omschakelingen. Helaas was dit model te complex om met de beschikbare rekencapaciteit te gebruiken. Dit model wordt beschreven in appendix A, meer duiding bij de werking van het model is terug te vinden in appendix B.
5.1 5.1.1
Modelformulering Indices Index
Beschrijving
j
Job: j ∈ J = {0, ..., N J}
k
Positie op machine: k ∈ K = {0, ..., N K}
m
Machine: m ∈ M = {0, ..., N M }
Tabel 5.1: Indices gebruikt in optimalisatiemodel
33
Hoofdstuk 5. Exact optimalisatiemodel
5.1.2
Variabelen Variable
Beschrijving
xjkm
= 1 als job j op positie k van machine m gepland is, anders 0.
lj
= 1 as job j te laat is, anders 0.
Tabel 5.2: Binaire variabelen gebruikt in optimalisatiemodel
Variable
Beschrijving
tj
Aantal dagen dat job j te laat is.
Tabel 5.3: Integer variabelen gebruikt in optimalisatiemodel
Variable
Beschrijving
cj
Tijdstip waarop job j wordt opgeleverd.
λj
Vaste kosten gemaakt door laattijdig leveren van job j.
Tabel 5.4: Continue variabelen gebruikt in optimalisatiemodel
34
Hoofdstuk 5. Exact optimalisatiemodel
5.1.3
35
Parameters Parameter
Beschrijving
Fjm
Vaste kost voor job j indien geproduceerd op machine m.
Vj
Variabele kost per dag dat job j te laat wordt opgeleverd.
Pjm
Productietijd in dagen van job j op machine m.
Rj
Tijdstip waarop job j in het systeem wordt toegevoegd.
Dj
De deadline van job j in dagen.
Sji
Tijdsduur in dagen nodig voor omschakeling van j naar i.
M
Een groot positief nummer.
NJ
Totaal aantal jobs dat ingepland moet worden.
NK
Maximaal aantal posities op een machine.
NM
Aantal machines.
Tabel 5.5: Parameters gebruikt in optimalisatiemodel
5.1.4
Doelfunctie M in.
NJ X
[(λj ) + (Vj ∗ tj )]
(5.1)
j=0
5.1.5
Randvoorwaarden N M X NK X
xjkm = 1,
j∈J
(5.2)
m=0 k=0 NJ X
xjkm ≤ 1,
k ∈ K, m ∈ M
(5.3)
j=0 NJ X j=0
xjkm −
NJ X
xj(k−1)m ≤ 0,
k ∈ K\{0}, m ∈ M
(5.4)
j=0
cj ≥ Pjm ∗ xjkm ,
j ∈ J, m ∈ M, k = 0
ci + M (2 − xikm − xj(k−1)m ) ≥ Pim + Sji + cj ,
j ∈ J, i ∈ J, k ∈ K\{0}, m ∈ M
(5.5)
(5.6)
Hoofdstuk 5. Exact optimalisatiemodel
cj ≥ Rj +
N M X NK X
36
(Pjm ∗ xjkm ),
j∈J
(5.7)
m=0 k=0
(cj − Dj ) ≤ M lj ,
j∈J
tj ≥ (cj − Dj ) − (1 − lj )M,
λj ≥ Fjm ∗ lj − M (1 −
NK X
j∈J
j ∈ J, m ∈ M
xjkm )
(5.8)
(5.9)
(5.10)
k=0
xjkm , lj ∈ {1, 0},
j ∈ J, i ∈ J, k ∈ K, m ∈ M
tj , cj , λj ≥ 0,
5.1.6
j∈J
(5.11)
(5.12)
Beknopte duiding bij vergelijkingen
5.1 De doelfunctie minimaliseert de som van de variabele en de vaste kosten die worden gemaakt wanneer een job te laat wordt opgeleverd. De parameter λ bevat de vaste kosten, die tevens rekening houden met de locatie waarop geproduceerd wordt. De variabele kosten zijn vervat in het product van de variabele kost per dag en het aantal dagen dat te laat wordt opgeleverd (= de variabele tj ). 5.2 Zorgt ervoor dat alle jobs die gepland moeten worden exact 1 positie op een van de machines worden toegewezen. 5.3 Op een bepaalde positie van een machine mag maximaal 1 job worden gepland. (Voorkomt dat er meerdere jobs op eenzelfde positie van eenzelfde machine ingepland worden.) 5.4 Er mogen geen gaten in de planning van een bepaalde machine zitten. Indien dit wel het geval is zorgt dit voor problemen voor het bepalen van de omschakelingen die plaatsvinden. Concreet betekent dit dat wanneer er een job gepland is op positie k van machine m, er tevens een job gepland moet zijn op positie k-1 van dezelfde machine. (Met uitzondering van de jobs die op de eerste positie van een machine ingepland staan.)
Hoofdstuk 5. Exact optimalisatiemodel
37
5.5 De eindtijd van de job die op de eerste positie gepland staat is gelijk aan zijn productietijd op die machine. 5.6 De eindtijd van de job is gelijk aan de eindtijd van de voorgaande job op dezelfde machine, plus de tijd nodig om de machine om te schakelen en de tijd nodig voor het produceren van de job zelf. 5.7 Een job kan niet worden geproduceerd alvorens deze in het systeem wordt opgenomen. (Release time constraint.) 5.8 Waardetoekenning van de lateness parameter, deze moet gelijk zijn aan 1 wanneer een job te laat wordt opgeleverd en gelijk aan 0 wanneer dit niet het geval is. 5.9 Waardetoekenning voor de tardiness parameter, die het aantal dagen dat een bepaalde job te laat geleverd wordt weergeeft. (Deze informatie is nodig voor het berekenen van de variabele kosten.) 5.10 Waardeberekening van de lambda variabele in de doelfunctie. De reden voor het introduceren van deze variabele is het lineariseren van deze doelfunctie: M in.
N PJ
[(
N M N K P P
[Fjm ∗ xjkm ] ∗ lj ) + (Vj ∗ tj )]
j=0 m=0 k=0
Het probleem met deze doelfunctie is dat het linkerdeel van de vergelijking een vermenigvuldiging van twee variabelen is en bijgevolg niet lineair is. Door het toevoegen van de extra variabele λ en de beperking 5.10 wordt dit probleem opgelost. 5.11 Geeft aan welke variabelen binair zijn. 5.12 Geeft aan welke variabelen positief (maar niet binair) moeten zijn.
5.2
Resultaten
In deze sectie wordt eerst kort gesproken over de parametersets die gebruikt werden om het model te testen, waarna de resulaten van de optimalisatie worden besproken. De resultaten die hier worden bekomen hebben vooral als nut de vergelijking met de resultaten van de heuristieken in het volgende hoofdstukken 7 en 8.
Hoofdstuk 5. Exact optimalisatiemodel
5.2.1
38
Breedte onderzoek
Parametersets gebruikt als input Een eerste groep parametersets is gegenereerd met als doel een zo breed mogelijke weergave van de verschillende problemen weer te geven. In dit geval werd niet uitgegaan van de realistische opstelling van het probleem, omdat er hier onmiddelijk sprake is van een groter aantal machines (32). Het aantal jobs varieert hier van 5 tot en met 25 (increment per 5 jobs), en het aantal machines varieert van 2 tot en met 10 (increment per 2 machines). Alle mogelijke combinaties tussen de bovenstaande aantallen worden gemaakt, wat resulteert in 25 parametersets. Het resultaat is een groep parametersets die zowel zeer eenvoudige als complexere problemen bevat. Dit betekent dat met behulp van deze groep parameters kan worden onderzocht welke de grootste problemen zijn die aan de hand van dit model kunnen worden opgelost. Resultaten van het breedte onderzoek Zoals hieboven reeds aangegeven is het belangrijkste doel hier een verkenning van de mogelijkheden van het exact optimalisatiemodel. Tijdens een eerste test werd ingesteld dat de optimalisatie stopt indien na 4000 seconden nog geen exacte oplossing was gevonden. De relevante resultaten worden weergegeven in grafiek 5.1. Uit deze grafiek kan worden afgeleid dat de meest interessante problemen voor het vergelijken van de performantie van heuristieken en exacte oplossingsmethodes zich qua grootteorde bevinden in het plannen van een 15 tot 20 tal jobs. Ook het effect van het aantal machines op de complexiteit van het probleem moet nogmaals worden onderzocht om hiervan een duidelijker beeld te krijgen.
5.2.2
Diepte onderzoek
Parametersets gebruikt als input Op basis van het breedte onderzoek werden eerste 50 parametersets gegenereerd die 15 jobs inplannen op 2 tot 10 machines (interval = 2 machines). De maximale berekeningstijd wordt vastgelegd op 8000 seconden. De resultaten van deze berekening worden grafisch voorgesteld in figuur 5.2
Hoofdstuk 5. Exact optimalisatiemodel
39
Sum of computationTime
4.000
# Seconden nodig voor oplossing
3.500
3.000
2.500
2.000
1.500
1.000
500
0 2 4 6 8 10
5 0 1
10 524 41
15 142 4.000
20 4.000 4.000
25 4.000 4.000
1 1
16 6
4.001 1.737
4.000 4.000
4.001 4.001
1
10
1.114
4.001
4.003
NJ
Figuur 5.1: Weergave van rekentijd nodig voor oplossen van verschillende job en machine combinaties met exact optimalisatiemodel, optimalisatie stopgezet indien geen resultaat na 4000 seconden.
Resultaten van het diepte onderzoek Uit deze resultaten blijkt dat slechts 24% van deze willekeurige problemen met 15 jobs opgelost konnen worden binnen de vooropgestelde tijdslimiet. Enkel wanneer er zes of meer machines waren om de jobs op te plannen bleek het mogelijk om een exacte oplossing te bekomen. Op basis van deze resultaten werd besloten om naast een beperkt aantal van deze complexe problemen ook een groter aantal eenvoudiger problemen op te nemen om te vergelijken met de heuristieke oplossingen. Dit houdt in dat eerst gezocht wordt naar de grens vanaf wanneer een probleem (gemiddeld genomen) te complex wordt om snel optimaal te kunnen oplossen. Om te onderzoeken hoe groot de problemen zijn die binnen een redelijk tijdsbestek kunnen worden opgelost worden 150 parametersets gegenereerd. Deze sets combineren 10 tot 15 jobs met 2 tot 10 machines. (Opnieuw worden machines per 2 toegevoegd.) Als conditie werd hier vooropgesteld dat de problemen binnen 1000 seconden tot een oplossing moeten komen. Indien deze tijd wordt overschreden wordt de oplossing voortijdig afgebroken. De resultaten van deze berekeningen worden weergegeven in figuur 5.3. Het is duidelijk dat
Hoofdstuk 5. Exact optimalisatiemodel
40
Sum of Solved
10
9
8
# Gevonden oplossingen
7
# Machines
6
2 4
5
6 8
4
10
3
2
1
0 # Jobs
Figuur 5.2: Exacte oplossingen gevonden binnen 8000 seconden voor het plannen van 15 jobs op M machines, per machinecategorie werden tien datasets getest.
aan de linkerzijde van de grafiek het grootste deel van de problemen binnen de 1000 seconden kan worden opgelost. Vanaf het moment dat er 12 (of meer) jobs ingepland moeten worden ontstaan er problemen indien er minder dan zes machines zijn. Een meer gedetailleerd beeld wordt gegeven in tabel 5.6. Het is duidelijk dat voor het plannen van 10 jobs de afwijkingen gemiddeld gezien zeer klein zijn, deze problemen zijn dus meestal binnen een zeer gering tijdsbestek op te lossen. Ook voor het plannen van een groter aantal jobs op 6 of meer machines zijn geen extreme evoluties merkbaar. De grote toename in complexiteit wordt echter gevonden wanneer meer dan 10 jobs op slechts twee machines moeten worden ingepland. Het is duidelijk dat de afwijkingen hier snel vrij groot blijven na 1000 seconden oplossingstijd.
Hoofdstuk 5. Exact optimalisatiemodel
41
Sum of Field1
100%
90%
Fractie opgelost binnen 1000 seconden
80%
70%
60% # Machines 2
50%
4 6
40%
8 10
30%
20%
10%
0% 10
11
12
13
14
15
# Jobs # Jobs
Figuur 5.3: Probleemgroottes die binnen een kort tijdsbestek exact kunnen worden opgelost.
# Jobs
# Machines 2
4
6
8
10
10
3%
19%
0%
0%
1%
11
129%
0%
22%
0%
0%
12
2820%
21%
11%
0%
0%
13
272%
125%
11%
1%
0%
14
2646%
24%
0%
16%
2%
15
100%
35%
46%
17%
7%
Tabel 5.6: Gemiddelde relatieve afwijking van de beste oplossing ten opzichte van de best gekende lower bound, voor verschillende combinaties van machines en jobs.
Hoofdstuk 6
Exacte oplossingsmethode: Illustraties en duiding Dit hoofstuk geeft meer uitleg bij de werking van het model uit hoofdstuk 5. Het is vooral bedoeld om mensen die minder ervaring hebben met de interpretatie van MIP optimalisatiemodellen zicht te geven op de redeneringen achter de vergelijkingen. Dit hoofdstuk beschrijft eerst de verschillende indices, variabelen en parameters die worden gebruikt in het model. Daarna worden alle randvoorwaarden grondig besproken en waar nodig met een voorbeeld ge¨ıllustreerd.
6.1
Indices
Deze sectie beschrijft de indices die gebruikt worden in het exacte optimalisatie model. Een belangrijke bemerking is hier dat in de werkelijke implementatie evenals de modelformulering in hoofdstuk 5 de indices steeds starten vanaf de waarde nul, waar in dit hoofdstuk de indices steeds vanaf de waarde 1 zullen starten. De reden hiervoor is dat de implementatie in C++ minder foutgevoelig is wanneer indices starten vanaf nul (de eerste positie op een array krijgt immers steeds als adres nul mee). Tellen vanaf nul is voor sommige mensen contra-intu¨ıtief, daarom is ervoor gekozen om in dit hoofdstuk alle indices vanaf de waarde 1 te laten starten. De j index wordt gebruikt om de verschillende jobs die ingepland moeten worden te identificeren. De definitie van een job is hier een bepaalde hoeveelheid van een bepaald product die aan een specifieke klant beloofd is voor een bepaalde datum. Deze index loopt van 1 tot aan het nummer van de laatste job. 42
Hoofdstuk 6. Exacte oplossingsmethode: Illustraties en duiding
43
De m index wordt gebruikt om een onderscheid te maken tussen de verschillende machines. Verschillende machines hebben verschillende productiemogelijkheden, verschillende productiesnelheden en een verschillende locatie. (Dit verschil in locatie is belangrijk voor het bepalen van de transportkost tot bij de klant.) De k index duidt de posities op de verschillende machines aan waarop een job ingepland kan worden. Om alle mogelijke productieplanningen toe te staan zal deze index lopen van 1 tot het aantal jobs dat in het model opgenomen is. (Op deze manier kan het meest extreme scenario waarin alle jobs op eenzelfde machine worden ingepland beschouwd worden.)
6.2
Beslissingsvariabelen
Hieronder worden de binaire, integer en continue variabelen uit het model kort besproken en waar nodig ge¨ıllustreerd.
6.2.1
De positie selectie variabele
De binaire xjkm variabele geeft aan of een bepaalde job op een bepaalde positie van een bepaalde machine wordt ingepland. Doordat er in dit geval sprake is van variabele omschakeltijden tussen producten is het eenvoudiger om een job een bepaalde positie toe te kennen dan de job een tijdstip mee te geven om te beginnen en/of te eindigen. Veronderstel een situatie waarin vier jobs geproduceerd moeten worden met twee machines. De planning van deze productie kan worden weergegeven met de xjkm variabelen uit tabel 6.1. De waarden van de variabelen de tabel 6.1 beschrijven een situatie waarin de eerste machine de eerste job produceert en de tweede machine de tweede, vierde en derde job produceert in die volgorde. Hier kan opgemerkt worden dat in elke “blok” van variabelen slechts een variabele de waarde 1 kan toegewezen worden. (Indien dit niet het geval is worden meerdere jobs op eenzelfde positie ingepland.) Op de manier waarop dit wordt afgedwongen wordt later teruggekeerd wanneer de randvoorwaarden die opgelegd worden aan het model worden besproken.
Hoofdstuk 6. Exacte oplossingsmethode: Illustraties en duiding
Positie 1
Positie 2
Positie 3
Positie 4
44
Machine 1
Waarde
Machine 2
Waarde
x111
1
x112
0
x211
0
x212
1
x311
0
x312
0
x411
0
x412
0
x121
0
x122
0
x221
0
x222
0
x321
0
x322
0
x421
0
x422
1
x131
0
x132
0
x231
0
x232
0
x331
0
x332
1
x431
0
x432
0
x141
0
x142
0
x241
0
x242
0
x341
0
x342
0
x441
0
x442
0
Tabel 6.1: Alle xjkm variabelen nodig om 4 jobs in te plannen op 3 machines, inclusief toegekende waarden.
Hoofdstuk 6. Exacte oplossingsmethode: Illustraties en duiding
6.2.2
45
De lateness variabele
De binaire variabele lj heeft de waarde 1 wanneer een job later dan de beloofde datum wordt opgeleverd. De beloofde datum die in het model wordt gebruikt is de beloofde datum min een forfait die wordt gerekend voor de productiestappen die na het breiproces nog moeten gebeuren. Daar het breiproces dat in beschouwing wordt genomen door het model de bottleneck van het proces is, en de productietijden van de andere stappen vrij voorspelbaar zijn, is dit een aanvaardbare assumptie.
6.2.3
Tardiness variabele
De tardiness variable tj is een integer variabele die het aantal dagen weergeeft die een bepaalde job j te laat wordt opgeleverd. Deze variabele wordt gelijk gesteld aan nul wanneer de job in kwestie op de beloofde tijd of vroeger wordt opgeleverd. Anders is deze variabele gelijk aan het verschil tussen de beloofde datum en de datum van afwerking, dat naar boven wordt afgerond.
6.2.4
Tijdstip van opleveren
De continue variabele cj geeft het tijdstip weer waarop de job wordt afgeleverd.
6.2.5
Vaste kosten voorgesteld door λ
De lambda is een hulpvariabele die in de doelfunctie de vaste kosten voorstelt. De reden waarom hiervoor een hulpvariabele gebruikt moet worden is dat de locatie van de productie afgeleid wordt uit de machine die voor de productie wordt gebruikt. Dit houdt in dat twee variabelen vermenigvuldigd worden in de onderstaande doelfunctie: M in.
N PJ
[(
N M N K P P
[Fjm ∗ xjkm ] ∗ lj ) + (Vj ∗ tj )]
j=1 m=1 k=1
Het probleem met deze doelfunctie is dat in het linkerdeel van de vergelijking de variabelen xjkm en lj worden vermenigvuldigd. Door het toevoegen van de extra variabele λ en randvoorwaarde 5.10 wordt dit probleem opgelost.
6.3
Parameters
Het gebruik van de verschillende parameters in het model wordt ge¨ıllustreerd aan de hand van een voorbeeld. Tabel 6.2 geeft de vaste kosten, variabele kosten, productietijden en deadlines
Hoofdstuk 6. Exacte oplossingsmethode: Illustraties en duiding
46
voor tien jobs weer. Hoe deze waarden worden berekend kan teruggevonden worden in 4.2. De tabel 6.2 toont dat afhankelijk van de machine die wordt gekozen om een bepaald product te produceren de vaste kosten en de productietijden verschillen. Wanneer de productietijd als waarde -1 heeft betekent dit dat de desbetreffende machine job j niet kan produceren. Bijgevolg zal een andere machine geselecteerd moeten worden voor deze job. In tabel 6.3 wordt voor de tien jobs de matrix met de omschakeltijden weergegeven. Dit zijn de waarden die worden voorgesteld door parameter Sji . Parameter S39 krijgt in deze tabel de waarde 8 toegekend. Aangezien de index t met tijdsintervallen van 15 minuten werkt betekent dit dat een technicus twee uur (= 8 * 15 min) nodig heeft om deze omschakeling te vervolledigen. Hoe deze parameters gebruikt worden in het model zal worden ge¨ıllustreerd in sectie 6.4.
6.4 6.4.1
Randvoorwaarden opgelegd aan het model Afdwingen van een werkbaar schema
randvoorwaarden 5.2, 5.3 en 5.4 zorgen ervoor dat de planning van de verschillende jobs correct gebeurt. Dit door af te dwingen dat elke job een positie krijgt op een machine, ervoor te zorgen dat op elke positie maximaal een job gepland staat en dat er geen gaten voorkomen in de planning. De werking van de verschillende vergelijkingen wordt hieronder ge¨ıllustreerd. randvoorwaarde 5.2 Vergelijking 5.2 zorgt ervoor dat elke job een positie op een machine krijgt toegewezen. Veronderstel dat er twee jobs gepland moeten worden op twee machines dan kunnen alle xjkm variabelen worden voorgesteld in tabel 6.4. Deze acht variabelen beschrijven alle mogelijke posities waarop deze twee jobs ingepland kunnen worden. De variabele x122 die de waarde 1 toegekend krijgt in deze tabel heeft aldus als betekenis dat de eerste job op de tweede positie van de tweede machine wordt geproduceerd. Wanneer vergelijking 5.2 wordt uitgeschreven voor de bovenstaande situatie wordt het volgende resultaat bekomen:
x111 + x112 + x121 + x122 = 0 + 0 + 0 + 1 = 1
(6.1)
Hoofdstuk 6. Exacte oplossingsmethode: Illustraties en duiding
47
j
Fj1
Fj2
Fj3
Vj
Pj1
Pj2
Pj3
Dj
1
81,23
947,48
947,48
2,70
2,29
1,97
-1,00
2
2
827,48
203,48
203,48
2,70
0,15
0,13
0,12
2
3
631,31
480,38
480,38
1079,69
0,99
0,85
0,79
2
4
721,13
546,98
546,98
1079,69
2,37
2,03
1,89
2
5
83,48
1007,48
1007,48
113,95
3,40
-1,00
-1,00
2
6
83,48
1007,48
1007,48
1079,69
7,40
6,34
5,92
2
7
514,52
1460,12
1460,12
1079,69
8,49
7,28
-1,00
2
8
258,08
175,79
175,79
1079,69
-1,00
-1,00
0,09
2
9
426,95
1195,25
1195,25
1079,69
0,08
0,07
-1,00
2
10
358,04
1151,48
1151,48
113,95
-1,00
-1,00
3,17
2
Tabel 6.2: Mogelijke invulling voor de waarden van de parameters in het optimalisatiemodel.
1
2
3
4
5
6
7
8
9
10
1
2
2
7
5
6
2
32
2
7
2
2
2
2
10
10
10
2
32
2
12
2
3
7
10
2
7
7
2
32
2
8
2
4
5
10
7
2
6
2
32
2
7
2
5
6
10
7
6
2
2
32
2
6
2
6
2
2
2
2
2
2
32
2
12
2
7
32
32
32
32
32
32
2
32
32
32
8
2
2
2
2
2
2
32
2
12
2
9
7
12
8
7
6
12
32
12
2
2
10
2
2
2
2
2
2
32
2
2
2
Tabel 6.3: Omschakeltijden in aantal tijdsintervallen van 15 minuten.
Hoofdstuk 6. Exacte oplossingsmethode: Illustraties en duiding
48
Machine 1
Waarde
Machine 2
Waarde
x111
0
x112
0
x211
0
x212
0
x121
0
x122
1
x221
0
x222
0
Positie 1
Positie 2
Tabel 6.4: Invloed van randvoorwaarde 5.2 op waarde van xjkm variabelen.
x211 + x212 + x221 + x222 = 0 + 0 + 0 + 0 6= 1
(6.2)
De berekeningen tonen dat vergelijking 6.1 gerespecteerd wordt, wat betekent dat job 1 aan exact 1 positie werd toegewezen. Vergelijking 6.2 daarentegen is in dit geval gelijk aan nul, wat betekent dat deze job op geen enkele machine ingepland staat. Dit is een situatie die in de realiteit niet werkbaar is, daar verondersteld wordt dat alle orders die in het model worden ingevoerd ook geproduceerd moeten worden. Het model zal er aldus voor zorgen dat de tweede job ook op een van de drie resterende mogelijke posities wordt ingepland. Wanneer variabele x111 ook de waarde 1 toegekend krijgt wordt vergelijking 6.1 gelijk aan 2. Net zoals een nulwaarde wordt dit eveneens niet aanvaard, het heeft immers geen nut om een order meermaals te produceren. randvoorwaarde 5.3 De waarden uit tabel 6.5 voldoen aan de conditie gesteld in randvoorwaarde 5.2. Dit schema is echter ook niet werkbaar in de praktijk, aangezien beide jobs als eerste van
Positie 1
Positie 2
Machine 1
Waarde
Machine 2
Waarde
x111
1
x112
0
x211
1
x212
0
x121
0
x122
0
x221
0
x222
0
Tabel 6.5: Invloed van randvoorwaarde 5.3 op waarde van xjkm variabelen.
Hoofdstuk 6. Exacte oplossingsmethode: Illustraties en duiding
49
start moeten gaan op de eerste machine. Deze machine kan echter niet simultaan twee jobs produceren dus moet het model verplicht worden om slechts 1 job per positie toe te wijzen. Dit gebeurt aan de hand van vergelijking 5.3 die hieronder wordt uitgeschreven voor dit voorbeeld:
k = 1, m = 1 :
x111 + x211 = 1 + 1 = 2 ≥ 1
(6.3)
k = 2, m = 1 :
x121 + x221 = 0 + 0 = 0 ≤ 1
(6.4)
k = 1, m = 2 :
x112 + x212 = 0 + 0 = 0 ≤ 1
(6.5)
k = 2, m = 2 :
x122 + x222 = 0 + 0 = 0 ≤ 1
(6.6)
De berekening 6.3 toont dat de randvoorwaarde die wordt opgelegd door vergelijking 5.3 hier niet wordt gerespecteerd. Wanneer het model ervoor zorgt dat deze randvoorwaarde gerespecteerd wordt zullen nooit twee jobs op een zelfde positie en machine ingepland zijn. randvoorwaarde 5.4 Om een correct schema af te dwingen moet tevens een meer technisch probleem worden voorkomen, de bovenstaande randvoorwaarden kunnen perfect vervuld zijn wanneer er gaten in de planning zitten zoals weergegeven in tabel 6.6.
Positie 1
Positie 2
Positie 3
Machine 1
Waarde
Machine 2
Waarde
x111
1
x112
0
x211
0
x212
0
x121
0
x122
0
x221
0
x222
0
x131
0
x132
0
x231
1
x232
0
Tabel 6.6: Invloed van randvoorwaarde 5.4 op waarde van xjkm variabelen.
Hoofdstuk 6. Exacte oplossingsmethode: Illustraties en duiding
50
Het probleem dat de situatie hierboven oplevert wordt veroorzaakt doordat er op de tweede positie van de eerste machine geen job ingepland staat, maar dit op de derde positie wel het geval is. Omdat de posities op de machines worden gebruikt voor het berekenen van de omschakeltijden tussen de verschillende jobs is dit schema binnen het model niet werkbaar. (Niettegenstaande dat in de praktijk deze lege positie gewoon zou kunnen worden overgeslaan en weinig problemen met zich mee zou brengen.) randvoorwaarde 5.4 zorgt ervoor dat bovenstaande situatie niet kan voorkomen. Hieronder wordt de vergelijking uitgeschreven voor de situatie in tabel 6.6 (enkel voor machine 1):
k = 2, m = 1 :
k = 3, m = 1 :
[x121 + x221 ] − [x111 + x211 ] = [0 + 0] − [1 + 0] = −1 ≤ 0
(6.7)
[x131 + x231 ] − [x121 + x221 ] = [0 + 1] − [0 + 0] = 1 ≥ 1
(6.8)
Het is duidelijk dat vergelijking 5.4 niet toestaat dat er gaten in de planning voorkomen (vergelijking 6.8 schendt immers de voorwaarde opgelegd door deze vergelijking). Een lege positie op een machine kan niet gevolgd worden door een positie waarop een job ingepland staat. Dit zorgt ervoor dat alle jobs gegroepeerd worden aan het begin van het schema. Doordat de jobs telkens op aansluitende posities gepland staan levert dit ook geen problemen op voor de berekening van de omschakeltijden.
6.4.2
Tijdstip bepalende randvoorwaarden
De randvoorwaarden die hier worden besproken zorgen ervoor dat het tijdstip waarop de productie van een job wordt afgerond correct wordt berekend. Het gaat het over vergelijkingen 5.5 tot en met 5.10. randvoorwaarde 5.5 randvoorwaarde 5.5 zorgt ervoor dat alle jobs die op de eerste positie van een machine gepland zijn een eindtijd meekrijgen die gelijk is aan de productietijd op de machine waarop ze gepland zijn. Het onderstaande voorbeeld verklaart de werking van deze vergelijking. In tabel 6.8 is opnieuw een situatie geschetst waar twee jobs op twee machines ingepland moeten worden. In dit geval staan deze jobs op de eerste positie van de eerste respectievelijk tweede machine. Dit betekent dat hun eindtijden bepaald worden aan de hand van vergelijking
Hoofdstuk 6. Exacte oplossingsmethode: Illustraties en duiding
Positie 1
Positie 2
51
Machine 1
Waarde
Machine 2
Waarde
x111
1
x112
0
x211
0
x212
1
x121
0
x122
0
x221
0
x222
0
Tabel 6.7: Input voor illustratie van randvoorwaarde 5.5
5.5. Wanneer deze vergelijkingen worden uitgeschreven en verondersteld wordt dat het hier om dezelfde job 1 en 2 gaat als wordt beschreven in tabel 6.2 wordt het volgende resultaat bekomen:
c1 ≥ P11 ∗ x111 ⇒ c1 ≥ 2.29 ∗ 1 ⇒ c1 ≥ 2.29
(6.9)
j = 1, m = 2,
c1 ≥ P12 ∗ x112 ⇒ c1 ≥ 1.97 ∗ 0 ⇒ c1 ≥ 0
(6.10)
j = 2, m = 1,
c2 ≥ P21 ∗ x211 ⇒ c2 ≥ 0.15 ∗ 0 ⇒ c2 ≥ 0
(6.11)
c2 ≥ P22 ∗ x212 ⇒ c2 ≥ 0.13 ∗ 1 ⇒ c2 ≥ 0.13
(6.12)
j = 1, m = 1,
j = 2, m = 2,
De berekeningen hierboven geven duidelijk aan dat de randvoorwaarde wordt opgelegd aan een job enkel en alleen als die job op de eerste positie van de machine staat (6.9 en 6.12). Indien dit niet het geval is (6.10 en 6.11) dwingt de randvoorwaarde enkel af dat het tijdstip waarop de job afgewerkt is groter is dan nul. randvoorwaarde 5.6 De randvoorwaarde 5.6 stelt dat het tijdstip waarop een job i afgerond is (ci ) groter moet zijn dan de som van het tijdstip waarop job j afgerond is (cj ), de tijd die nodig is om de omschakeling van j naar i uit te voeren (sji ) en de tijd nodig om job i te produceren op machine m (Pim ) als en slechts als job j net voor job i gepland is op eenzelfde machine m. Deze randvoorwaarde wordt hier opnieuw ge¨ıllustreerd aan de hand van een eenvoudig voorbeeld in tabel 6.8. Opnieuw wordt de situatie verondersteld waarin twee jobs gepland moeten
Hoofdstuk 6. Exacte oplossingsmethode: Illustraties en duiding
52
worden op twee machines. De situatie die hieronder beschreven wordt voldoet aan alle restricties beschreven in 6.4.1. De jobs die gepland moeten worden gespecificeerd als jobs 1 en 2 uit tabel 6.2 en tabel 6.3. Vergelijking 5.6 kan als volgt worden uitgeschreven (enkel machine 1 wordt in beschouwing genomen):
i = 1, j = 1, k = 2, m = 1,
c1 + M ∗ (2 − x121 − x111 ) ≥ P11 + s11 + c1
c1 + M ∗ (2 − 0 − 1) ≥ P11 + s11 + c1 c1 + M ≥ P11 + s11 + c1 (6.13)
i = 2, j = 1, k = 2, m = 1,
c2 + M ∗ (2 − x221 − x111 ) ≥ P21 + s12 + c1
c2 + M ∗ (2 − 1 − 1) ≥ P21 + s12 + c1 c2 + 0 ≥ 1.97 + 0.02 + 2.29 (6.14)
Vergelijking 6.13 heeft geen invloed op de berekening van het eind tijdstip. Dit doordat het hier niet gaat over twee jobs die na elkaar op eenzelfde machine ingepland worden. Dit wordt gecontroleerd door de uitdrukking M ∗ (2 − x121 − x111 ) die ervoor zorgt dat de big M waarde enkel wegvalt uit de vergelijking wanneer het gaat over twee jobs die na elkaar op eenzelfde machine ingepland staan.
Positie 1
Positie 2
Machine 1
Waarde
Machine 2
Waarde
x111
1
x112
0
x211
0
x212
0
x121
0
x122
0
x221
1
x222
0
Tabel 6.8: Input voor illustratie van randvoorwaarde 5.6
Hoofdstuk 6. Exacte oplossingsmethode: Illustraties en duiding
53
De tweede uitdrukking 6.14 gaat doordat de big M parameter wegvalt echter wel afdwingen dat het tijdstip waarop job nummer 2 afgerond is minimaal groter is dan het tijdstip waarop job 1 afgerond is (deze job was op de positie voor job 2 gepland op dezelfde machine), plus de tijd nodig voor de omschakeling, plus de tijd nodig voor het produceren van job 2 zelf. Zoals eerder vermeld wordt het tijdstip waarop de jobs op de eerste positie klaar zijn bepaald door vergelijking 5.5. Alle andere posities worden bepaald door deze vergelijking. randvoorwaarde 5.7 Deze randvoorwaarde geeft orders een release time mee. Dit betekent dat een order niet mag worden uitgevoerd voordat de release time van een order is aangebroken. Voor het optimalisatiemodel betekent dit dat het tijdstip waarop de job wordt afgeleverd (cj ) ook afhankelijk is van het tijdstip waarop een order in het systeem wordt opgenomen. Veronderstel dat order 1 een als release time R1 de waarde 5 meegekregen heeft, dit betekent dat de productie van order 1 niet mag starten voor dag 5. In de onderstaande vergelijking wordt verondersteld dat job 1 wordt geproduceerd op machine 2.
c1 ≥ R1 +
N M X NK X
(P1m ∗ x1km ) ⇒ c1 ≥ 5 + 1.97
(6.15)
m=1 k=1
De bovenstaande vergelijking zorgt ervoor dat het tijdstip waarop job 1 ten vroegste afgerond kan worden gelijk wordt aan 6.97. (Merk op dat het getal waarmee de release time wordt verhoogd de productietijd is op de machine waarop de desbetreffende job wordt geproduceerd.) Het is onnodig om bij de bovenstaande vergelijkingen de tijd nodig voor het omschakelen te betrekken, indien de omschakeling een invloed heeft op het starttijdstip wordt dit immers door randvoorwaarde (5.6) afgedwongen. randvoorwaarde 5.8 Deze randvoorwaarde zorgt ervoor dat de binaire variabele lj die aangeeft of een job al dan niet te laat zal worden opgeleverd de juiste waarde meekrijgt. Dit betekent dat wanneer de datum van oplevering cj later valt dan de beloofde datum Dj voor het desbetreffende order, de waarde van lj gelijkgesteld wordt aan 1.
Hoofdstuk 6. Exacte oplossingsmethode: Illustraties en duiding
54
Dit kan worden ge¨ıllustreerd door verder te bouwen op het voorbeeld van randvoorwaarde 5.6. Veronderstel dat de effectieve omschakeltijd sji in vergelijking 6.14 gelijk is aan 0.05. Dit zou dan betekenen dat het tijdstip waarop job 2 wordt opgeleverd gelijk zou zijn aan 4.31. (c2 = 4.31) Indien deze job beloofd was op dag 5 (D2 = 5) krijgen we de volgende invulling van vergelijking 5.8.
(c2 − D2 ) ≤ M ∗ l2 ⇒ 4.31 − 5 ≤ M ∗ lj
(6.16)
In de bovenstaande vergelijking kan de waarde voor lj zonder problemen op nul worden gezet. Omdat lj in de doelfunctie een bepaalde kost heeft meegekregen zal de optimalisatie van het model er tevens voor zorgen dat deze variabele ook op nul gezet zal worden indien dit mogelijk is. Wanneer dit order echter op dag 4 beloofd was krijgen we de onderstaande situatie:
(c2 − D2 ) ≤ M ∗ l2 ⇒ 4.31 − 4 ≤ M ∗ lj
(6.17)
Deze vergelijking zal in tegenstelling tot vergelijking 6.17 wel afdwingen dat lj de waarde 1 krijgt, daar het order later dan de beloofde datum wordt opgeleverd. randvoorwaarde 5.9 Omdat niet enkel het al dan niet te laat opgeleverd worden van een job van belang is voor de kosten toegekend aan een order moet er tevens nog een integer variabele worden gespecificeerd die aangeeft hoe veel dagen de order te laat werd opgeleverd. Dit gebeurt met behulp van de integer variabele tj die zijn waarde toegekend krijgt door deze randvoorwaarde. Hieronder een beknopt voorbeeld om de werking van deze vergelijking te illustreren. Een order j wordt opgeleverd op het moment 7.35, wanneer verondersteld wordt dat dit order beloofd was op dag 10 krijgen we de volgende invulling van deze randvoorwaarde:
tj ≥ (cj − Dj ) − (1 − lj ) ∗ M ⇒ tj ≥ (7.35 − 10) − (1 − 0) ∗ M ⇒ tj ≥ −M
(6.18)
De bovenstaande vergelijking laat op analoge wijze aan de vergelijkingen die de waarde voor de binaire lj variabele bepalen de mogelijkheid aan de variabele tj om gelijk aan nul te zijn. De waarde van tj heeft eveneens een impact op de doelfunctie, waardoor deze - indien mogelijk
Hoofdstuk 6. Exacte oplossingsmethode: Illustraties en duiding
55
- ook steeds de waarde nul zal meekrijgen. Wanneer deze job echter beloofd zou zijn op dag 4 krijgt de vergelijking de onderstaande vorm:
tj ≥ (cj − Dj ) − (1 − lj ) ∗ M ⇒ tj ≥ (7.35 − 4) − (1 − 1) ∗ M ⇒ tj ≥ 3.35
(6.19)
Hier zal de waarde van tj gelijk worden gesteld aan 4, omdat tj als een integer werd gedefinieerd. randvoorwaarde 5.10 Deze randvoorwaarde is opgenomen in het model om de correcte waarde voor de totale vaste kosten voor een bepaalde job j te gaan berekenen. De doelfunctie in zijn meest eenvoudige logische vorm kan als volgt worden geformuleerd: NJ N M X NK X X [Fjm ∗ xjkm ] ∗ lj ) + (Vj ∗ tj )] M in. [(
(6.20)
j=1 m=1 k=1
Niettegenstaande deze schrijfwijze conceptueel eenvoudiger is brengt deze problemen met zich K N M N P P [Fjm ∗ mee omdat het deel van deze functie dat de vaste kosten beschrijft, namelijk: ( m=1 k=1
xjkm ] ∗ lj ) niet lineair is. (De xjkm en lj variabelen worden met elkaar vermenigvuldigd.) Om deze vergelijking te lineariseren wordt de variabele λ ingevoerd die het totaal van de variabele kosten gaat voorstellen in de doelfunctie. Vergelijking 5.10 gaat aan λ de vaste kost toekennen van de machine waarop de job j wordt geproduceerd. Dit kan ge¨ıllustreerd worden aan de hand van een eenvoudig voorbeeld. Stel dat job 1 werd ingepland op machine 1 en te laat wordt opgeleverd, job 2 werd ingepland op machine 1 en wordt op tijd opgeleverd. (De waarden van de parameters worden opnieuw uit tabel 6.2 gehaald.
λ1 ≥ F11 ∗ l1 − M ∗ (1 −
NK X
x1k1 ) ⇒ λ1 ≥ 81.23 ∗ 1 − M ∗ (1 − 1) ⇒ λ1 ≥ 81.23
(6.21)
k=1
λ2 ≥ F21 ∗ l2 − M ∗ (1 −
NK X
x2k1 ) ⇒ λ2 ≥ 947.48 ∗ 0 − M ∗ (1 − 1) ⇒ λ2 ≥ 0
(6.22)
x1k1 ) ⇒ λ1 ≥ 81.23 ∗ 1 − M ∗ (1 − 0) ⇒ λ1 ≥ −M
(6.23)
k=1
λ1 ≥ F12 ∗ l1 − M ∗ (1 −
NK X k=1
Hoofdstuk 6. Exacte oplossingsmethode: Illustraties en duiding
56
Vergelijking 6.21 toont dat wanneer een job te laat wordt opgeleverd op een bepaalde machine, de job parameter λ de correcte waarde meekrijgt, zijnde de vaste kost geassocieerd met te laat produceren op de desbetreffende machine. (Noot: de reden waarom dit verschilt, is dat de locatie van de machine een impact heeft op de verzendingskosten van het product naar de klant die in dit geval gedragen moeten worden door het bedrijf zelf.) Vergelijking 6.23 toont dat wanneer de vergelijking wordt uitgeschreven voor een machine die de job niet produceert, de randvoorwaarde geen kost zal opleggen. Er wordt eveneens geen randvoorwaarde opgelegd wanneer een job op tijd wordt geproduceerd, dit wordt ge¨ıllustreerd in vergelijking 6.22.
6.4.3
Eigenschappen van variabelen bepalen
randvoorwaarden 5.11 en 5.12 bepalen de limieten en de aard van de verschillende beslissingsvariabelen. Deze vergelijkingen zijn zeer eenvoudig en worden hier niet verder besproken.
Hoofdstuk 7
Genetische metaheuristiek Een frequent gebruikte techniek voor het optimaliseren van planningen voor parallelle machines is een genetisch algoritme. Dit type algoritme kan voor een enorme verscheidenheid aan problemen worden aangewend. Hier wordt een genetisch algoritme toegepast op het probleem zoals het ook in hoofdstuk 5 werd behandeld. De metaheuristiek die hier wordt beschreven is in grote lijnen gebaseerd op Vallada & Ruiz (2011) en technieken besproken door Kurz & Askin (2001) en Ruiz & Allahverdi (2007). Genetische heuristieken zijn gebaseerd op biologische principes uit de evolutieleer van Darwin. Deze heuristieken maken gebruik van populaties waarin de best aangepaste leden meer kans hebben om zich voort te planten dan minder aangepaste leden. Deze nakomelingen erven eigenschappen van twee van hun ouders en hebben zo de kans om positieve eigenschappen van beide te gaan combineren. Na verloop van tijd gaan deze nieuwe leden de bestaande leden van de populatie gaan verdringen en gaan deze zelf over tot het produceren van nakomelingen in de vorm van gecombineerde oplossingen. Dit hoofdstuk bespreekt eerst de verschillende elementen van de gebruikte heuristiek. Het tweede deel van dit hoofdstuk analyseert de kalibratie van de verschillende parameters van dit model. In het derde en laatste deel van dit hoofdstuk worden de resultaten voor de verschillende probleemgroottes besproken.
57
Hoofdstuk 7. Genetische metaheuristiek
7.1
58
Werking van de metaheuristiek
Deze paragraaf geeft een uitgebreide uitleg bij de functionaliteiten en parameters waaruit deze metaheuristiek bestaat. Deze worden waar nodig ge¨ıllustreerd om de concepten zo duidelijk mogelijk voor te stellen.
7.1.1
Representatie van productieschemas
Een productieschema wordt weergeven in een multidimensionele array waarbij de eerste dimensie de machine voorstelt en de tweede dimensie de positie op de machine voorstelt. Deze representatie wordt ge¨ıllustreerd in figuur 7.1. Op deze figuur wordt onder andere aangeduid dat de negende job op de derde positie van de tweede machine ingepland staat.
Figuur 7.1: Representatie van genetisch algoritme
7.1.2
Creatie van de populatie
Voor het opstellen van een zo divers mogelijke populatie kunnen zowel volledig willekeurige technieken als enkelvoudige heuristieken worden gebruikt. Voorlopig wordt enkel de willekeurige generatie besproken, maar in sectie 7.3.1 wordt dieper ingegaan op enkele heuristieken die kwalitatief betere oplossingen toevoegen aan de initi¨ele populatie. Het is belangrijk om hierbij op te merken dat bij de generatie van deze schema’s wordt voorkomen dat niet werkbare schema’s worden geproduceerd. (Dit is immers mogelijk doordat niet alle producten door alle machines kunnen worden geproduceerd.) Dit is echter geen vereiste, een algoritme zou ook perfect kunnen starten vanuit een aantal niet volledig werkbare oplossingen. Het algoritme dat hiervoor wordt gebruikt bestaat uit de volgende stappen: 1. Selecteer een willekeurige job uit de nog niet ingeplande jobs. 2. Voor deze job wordt willekeurig een machine geselecteerd, een machine kan enkel worden geselecteerd wanneer deze de job kan produceren.
Hoofdstuk 7. Genetische metaheuristiek
59
3. Plan deze job in op de volgende vrije positie op deze machine.
7.1.3
Selectie van ouders
De techniek die wordt gebruikt voor het selecteren van de leden uit de populatie die nakomelingen is een tornooi selectie mechanisme. Dit houdt in dat een willekeurig geselecteerde groep van vooraf bepaalde omvang uit de populatie tegen elkaar wordt uitgespeeld in een tornooi, de winnaar van dit tornooi kan nakomelingen produceren. De groep die tegen elkaar wordt uitgespeeld is willekeurig bepaald met behulp van een parameter die pressure% genoemd wordt. Het percentage dat deze parameter voorstelt is de fractie van de populatie die meedoet aan het tornooi. De deelnemers worden telkens opnieuw willekeurig geselecteerd tot het percentage bereikt is. Voor het selecteren van twee ouders worden twee aparte tornooien gebruikt, de respectievelijke winnaars worden dan gekruist om nakomelingen te produceren. Uit de literatuur blijkt dat deze manier van selecteren vaak betere resultaten heeft dan de klassieke roulette selectie methode (Glover & Kochenberger, 2003).
7.1.4
Generatie van nakomelingen
De techniek die wordt gebruikt voor de generatie van nakomelingen is dezelfde als in Vallada & Ruiz (2011), namelijk een one-point crossover. Deze techniek wordt hier ge¨ıllustreerd met behulp van figuren 7.2 en 7.3.
Figuur 7.2: Ouders geselecteerd uit populatie met behulp van tournament selectie mechanisme.
Hoofdstuk 7. Genetische metaheuristiek
60
Figuur 7.3: Nakomelingen gengenereerd uit ouders met willekeurige crossover.
Figuur 7.2 geeft een voorbeeld van twee ouders in een situatie waarbij 13 jobs ingepland moeten worden op twee machines. Deze twee ouders zijn geselecteerd uit de populatie met behulp van de tournament selectie mechanismen beschreven in 7.1.3. Op elke machine van de eerste ouder wordt een willekeurige positie p gekozen gekozen uit de posities die gebruikt worden door jobs. Hier is dit de vierde positie op de eerste machine en de derde positie op de tweede machine. Vervolgens worden posities nul tot en met p gekopieerd naar de eerste nakomeling en posities p + 1 tot en met de laatste ingenomen positie worden gekopieerd naar de tweede nakomeling. Dit zijn de jobs die in het grijs worden weergegeven op figuur 7.3. De jobs die op dit moment nog niet in de nakomelingen voorkomen worden gepland in dezelfde volgorde en op dezelfde machine als ze in de tweede ouder voorkomen. Voor beide nakomelingen worden dus op alle machines van de tweede ouder alle posities sequentieel overlopen en indien de job nog niet gepland is op de machine wordt de job ingepland in de volgorde zoals die op de tweede ouder wordt waargenomen. Voor de eerste nakomeling zien we bijvoorbeeld dat job nummer 13 op de eerste positie van de tweede ouder nog niet is ingepland. Hierdoor krijgt deze job de eerstvolgende positie op deze machine. De jobs die worden toegevoegd aan de hand van informatie uit de tweede ouder worden weergegeven in het groen op figuur 7.3.
7.1.5
Mutaties
Niettegenstaande dat de omvang van de populate een zekere diversiteit garandeert in het begin van het algoritme is het toch mogelijk dat een genetisch algoritme naar een lokaal
Hoofdstuk 7. Genetische metaheuristiek
61
optimum evolueert. Om dit te voorkomen wordt het concept mutaties toegevoegd aan de crossover. Volgens een bepaalde waarschijnlijkheid (die vooraf wordt bepaald) vindt er een mutatie plaats in de nakomelingen, waardoor deze eigenschappen krijgen die niet noodzakelijk aanwezig waren bij hun ouders.
Figuur 7.4: Mutaties bij nakomelingen bevorderen de diversiteit in de populatie.
Dit principe kan natuurlijk op verschillende manieren worden ingevuld. Hier wordt zeer eenvoudig een volledig willekeurige job geselecteerd, die wordt verwijderd vanop de positie waar deze zich nu bevind en op een andere willekeurige positie opnieuw wordt ingepland. Om de berekeningen hier zo snel mogelijk te laten verlopen wordt er geen rekening gehouden met het het al dan niet mogelijk zijn van deze nieuwe volgorde. (Indien deze combinatie niet mogelijk is wordt deze sowieso afgestraft bij de generatieovergang omdat de doelfunctie nooit in de buurt kan komen van een werkbare oplossing.) Een interesante techniek om mutaties effectiever te gebruiken wordt voorgesteld door Sewell et al. (2008). Er wordt gesteld dat door leden van de populatie die hoger gerangschikt zijn een lagere kans te geven op mutaties bij hun nakomelingen en vice versa de resultaten van de genetische heuristiek kunnen worden verbeterd.
7.1.6
Generatieovergang
Wanneer het vooropgestelde aantal nakomelingen gegenereerd is begint een fase die generatieovergang genoemd kan worden. Tijdens deze fase wordt een deel van de huidige populatie vervangen door leden van de nieuwe generatie. Dit kan aan de hand van verschillende methoden, maar om zo veel mogelijk kwalitatieve oplossingen te bewaren is hier gekozen voor een methode die enkel nieuwe leden toelaat in de populatie als ze beter zijn dan het slechtste lid
Hoofdstuk 7. Genetische metaheuristiek
62
van de huidige populatie. Beter zijn betekent hier natuurlijk een betere (kleinere) waarde hebben voor de doelfunctie. Concreet worden alle nakomelingen sequenteel overlopen, indien de doelfunctie een betere waarde heeft dan die van het slechtste lid van de populatie wordt de populatie eerst doorzocht om te controleren of er geen identieke oplossing al is opgenomen in de populatie. (Indien deze controle niet gebeurt kan de diversiteit in de populatie zeer snel afnemen.) Wanneer geen identiek lid wordt gevonden in de populatie wordt het laatste lid van de populatie verwijderd en de oplossing ingevoegd op de correcte positie. Het vergelijken van het nieuwe slechtste lid van de populatie met de rest van de nakomelingen gaat dan op deze manier verder tot alle nakomelingen overlopen zijn en de nieuwe populatie volledig is. Deze nieuwe populatie kan dan opnieuw worden gebruikt om nakomelingen te genereren en zo dichter bij een optimale oplossing te komen.
7.2
Kalibratie van parameters
Om te garanderen dat een metaheuristiek optimaal presteert is het kalibreren van de parameters een van de belangrijkste stappen. In deze sectie worden de resultaten ten opzichte van de optimale oplossing vergeleken voor verschillende probleemgroottes en voor verschillende waarden van de parameters. De parameters die worden besproken zijn: het aantal leden in de populatie, het pressure% (dit is de fractie van de populatie die wordt geselecteerd voor de tornooi selectiemethode), de kans op mutaties, het aantal generaties en het aantal nakomelingen per generatie. Als input voor deze kalibratie werden drie verschillende probleemgroottes geselecteerd, een klein probleem dat 5 jobs inplant op twee machines, een middelgroot probleem dat 10 jobs inplant op 4 machines en een groot probleem dat 15 jobs inplant op 10 machines. Let wel dat deze problemen allemaal nog steeds klein zijn, om te garanderen dat deze nog steeds optimaal kunnen worden opgelost met behulp van het vereenvoudigde exact optimalisatiemodel. Op deze manier kan met zekerheid worden vastgesteld wanneer de optimale oplossing bereikt is. Aan het eind van deze sectie wordt tevens onderzocht hoe de parameters aangepast moeten worden om grotere problemen efficient en effectief te kunnen oplossen. Daar kan niet meer
Hoofdstuk 7. Genetische metaheuristiek
Populatie omvang
63
50
Pressure %
20%
Kans op mutatie
10%
Aantal generaties Aantal nakomelingen
10 100
Tabel 7.1: Basiswaarden voor parameters gebruikt tijdens kalibratie.
worden uitgegaan van de optimale oplossing omdat het probleem te groot is geworden om het met exacte technieken op te lossen. Voor het kalibreren worden alle parameters - met uitzondering van de parameter die gekalibreerd wordt - vastgezet op de waarden uit tabel 7.1. De resultaten van deze kalibratie worden steeds weergegeven in grafieken zoals figuur 7.5. De horizontale as van deze grafiek toont telkens de verschillende waarden voor de parameter die wordt gekalibreerd. De rechtse verticale as toont het gemiddelde eindresultaat in procent ten opzichte van de optimale waarde en is gelinkt aan de lijnen in de grafiek. De bedoeling is hier om zo dicht mogelijk bij de optimale waarde van 100% te komen. (Afwijkingen ten opzichte van deze waarde kunnen enkel naar boven daar dit een minimalisatieprobleem betreft.) De optimale waarden werden bekomen met behulp van het model uit hoofdstuk 5. De linkse verticale as geeft de rekentijd weer en is gelinkt met de staafgrafieken. Deze tijd wordt uitgedrukt in seconden en het is natuurlijk de bedoeling om deze rekentijd zo laag mogelijk te houden. De verschillende kleurschakeringen van de kolommen en lijnen geven dan weer de omvang van het probleem weer, er wordt steeds een klein, middengroot en groot probleem getest. (Al is het grootste probleem hier nog steeds vrij bescheiden qua omvang omdat de rekencapaciteit niet onbeperkt is en anders zeer lange rekentijden nodig zouden zijn voor deze kalibratie.) Na deze eerste tests volgt nog een test op grotere schaal die wordt beschreven in sectie 7.2.6.
Hoofdstuk 7. Genetische metaheuristiek
7.2.1
64
Omvang van de populatie
Op grafiek 7.5 worden de resultaten van de kalibratierun voor de populatieomvang weergegeven. Het is duidelijk dat voor het kleinste probleem een relatief kleine omvang van de populatie
180
180%
160
160%
140
140%
120
120%
100
100%
80
80%
60
60%
40
40%
20
20%
Objective function value relative to optimal value
Computation time in seconds
(50) lijkt te volstaan om het optimum met een vrij grote zekerheid te kunnen bereiken.
0%
0 20
40
60
80
100
120
140
160
180
Population size Time Small
Time Medium
Time Big
% of optimum Small
% of optimum Medium
% of optimum Big
Figuur 7.5: Kalibratie van populatieomvang
Voor het middelgroot probleem lijkt toch een iets grotere populatie van om en bij de 100 leden aangewezen. Voor het grootste probleem dat hier in beschouwing wordt genomen is er duidelijk zelfs bij grote omvang van de populatie nog steeds geen garantie dat het optimum wordt gevonden, maar vanaf om en bij de 100 leden in de populatie blijft de metaheuristiek hier steeds binnen 20% van de optimale waarde van de doelfunctie. Het is duidelijk dat de rekentijd slechts een geringe toename vertoont bij toenemende groote van de populatie. Dit gecombineerd met de waarneming dat een grotere populatie voor een beter resultaat kan gaan zorgen doet besluiten dat de populatiegrootte best op 100 of meer leden wordt vastgelegd. Voor grotere problemen zal het waarschijnlijk aan te raden zijn om deze parameter nog verder te verhogen om de diversiteit van de mogelijke oplossingen te kunnen garanderen. Deze toename zal slechts een geringe verhoging van de rekentijd veroorzaken, maar de resultaten
Hoofdstuk 7. Genetische metaheuristiek
65
naar alle waarschijnlijkheid positief be¨ınvloeden.
7.2.2
Pressure%
Grafiek 7.6 toont de resultaten voor de calibratie van de pressure% parameter die de competitie binnen de populatie om in aanmerking te komen voor reproductie aanduidt. De oplossingstijd die de metaheuristiek nodig heeft wordt blijkbaar niet be¨ınvloed door deze parameter, daar geen relevante verschillen of trends op te merken zijn in deze waarden. Met de rekentijd hoeft aldus geen rekening worden gehouden voor het kiezen van een waarde voor deze parameter. 120
250%
Computation time in seconds
200%
80 150%
60
100% 40
50%
Objective function value relative to optimal value
100
20
0
0% 5%
10%
15%
20%
25%
30%
35%
40%
45%
50%
Pressure% Time Small
Time Medium
Time Big
% of optimum Small
% of optimum Medium
% of optimum Big
Figuur 7.6: Kalibratie van pressure%
De eindoplossing van het algoritme wordt echter wel be¨ınvloed door de waarde van deze parameter. Wanneer de waarde van pressure% boven de 25% uitstijgt worden snel grote afwijkingen (> 20%) ten opzichte van het optimum waargenomen. Daarom lijkt het aangewezen om deze waarde binnen het interval tussen 10% en 20% te houden. Een hogere waarde van deze parameter zal betere oplossingen bevoordelen, maar zal ervoor zorgen dat een deel van de oplossingen nooit aan bod komt. (Het percentage van de oplossingen dat nooit aan bod komt is exact gelijk aan de waarde die gekozen wordt voor deze
Hoofdstuk 7. Genetische metaheuristiek
66
parameter). Aan het andere eind van het spectrum zal een waarde die te laag is goede oplossingen onvoldoende gaan bevoordelen wanneer deze moeten voortplanten. Dit heeft dan tot resultaat dat de nakomelingen van minder goede kwaliteit zijn.
7.2.3
Kans op mutaties
De resultaten voor de mutatiekans worden weergegeven in grafiek 7.7. Opnieuw is er geen relelevante variatie waar te nemen in de rekentijden voor verschillende waarden van deze parameter. 120
160%
140%
Computation time in seconds
120%
80 100%
60
80%
60% 40
40%
Objective function value relative to optimal value
100
20 20%
0
0% 5%
10%
15%
20%
25%
30%
35%
40%
45%
50%
Mutation probability in offspring Time Small
Time Medium
Time Big
% of optimum Small
% of optimum Medium
% of optimum Big
Figuur 7.7: Kalibratie van kans op mutatie bij nakomelingen
De invloed op de effectiviteit van het algoritme van deze parameter is tevens onduidelijk. De beste zone lijkt om en bij de 40% te zijn, onafhankelijk van probleemgrootte. Het gebruik van mutaties in genetische algoritmes wordt ook in de literatuur soms afgeschreven als zijnde van weinig nut, daar het selectiemechanisme de voornaamste werkende component is. (Glover & Kochenberger, 2003). Ondanks het geringe effectiviteit van de mutaties in dit voorbeeld wordt geopteerd om deze hier toch binnen het algoritme te houden. Deze keuze wordt gemaakt omdat de meerderheid
Hoofdstuk 7. Genetische metaheuristiek
67
van de auteurs deze stap toch nuttig achten. Tevens zorgt deze bewerking niet voor een merkbare vertraging van het algoritme en is er dus weinig reden om deze niet op te nemen.
7.2.4
Aantal generaties
Figuur 7.8 geeft de resultaten weer voor de kalibratie van het aantal generaties. Een eerste belangrijke conclusie is dat de rekentijd sterk toeneemt wanneer meer generaties worden gebruikt in de metaheuristiek. Deze relatie blijft echter lineair en zal veel minder snel dan de exacte optimalisatie gaan toenemen wanneer het probleem complexer wordt. De afwijking ten opzichte van de optimale oplossing wordt gemiddeld genomen kleiner wanneer er meer generaties worden gebruikt. Het is echter moeilijk om te bepalen wanneer precies het optimale aantal generaties is bereikt. Het bepalen van deze parameter wordt nog moeilijker voor grotere problemen. 140%
900
800
700
Computation time in seconds
100% 600
80%
500
400
60%
300 40% 200
Objective function value relative to optimal value
120%
20% 100
0
0% 10
20
30
40
50
60
70
80
90
# Generations Time Small
Time Medium
Time Big
% of optimum Small
% of optimum Medium
% of optimum Big
Figuur 7.8: Kalibratie van aantal generaties
Een mogelijkheid om deze problematiek op te lossen is het bijhouden van de verbetering in de opeenvolgende generaties. Wanneer bijvoorbeeld gedurende n generaties geen nieuwe nakomelingen meer werden opgenomen in de populatie kunnen de iteraties worden gestopt. Daarnaast moet natuurlijk ook een maximaal aantal iteraties of tijdslimiet worden vastgelegd. Hierop wordt later teruggekomen in sectie 7.3.2.
Hoofdstuk 7. Genetische metaheuristiek
68
Op basis van de resultaten in figuur 7.8 wordt gesteld dat meer dan 50 generaties naar alle waarschijnlijkheid geen drastische verbetering zal aanbrengen in de kwaliteit van de oplossing. (Of althans niet bij problemen van deze grootteorde.)
7.2.5
Aantal nakomelingen
Zoals grafiek 7.9 aantoont is het aantal nakomelingen per generatie opnieuw een parameter die een lineaire toename van de oplossingstijd veroorzaakt. Deze toename is echter veel minder steil dan de toename bij het aantal generaties. 300%
200
180
Computation time in seconds
140 200% 120
150%
100
80 100% 60
40
Objective function value relative to optimal value
250% 160
50% 20
0
0% 10
20
30
40
50
60
70
80
90
100
110
120
130
140
Offspring per generation Time Small
Time Medium
Time Big
% of optimum Small
% of optimum Medium
% of optimum Big
Figuur 7.9: Kalibratie van aantal nakomelingen
De waarde van deze variabele lijkt echter geen eenduidige relatie te hebben met het resultaat van de metaheuristiek. Het lijkt echter veilig om het aantal nakomelingen minimaal gelijk te stellen aan de omvang van de populatie. De impact van deze parameter wordt verder bestudeerd in sectie 7.2.6.
7.2.6
Toepassing op probleem met realistische grootte
Aanpassingen van de parameterwaarden voor grotere problemen Na de bovenstaande testen wordt de heuristiek uitgetest op een probleem van realistische omvang, 300 jobs inplannen op 32 machines. Hierbij wordt ook opnieuw onderzocht hoe
Hoofdstuk 7. Genetische metaheuristiek
69
gevoelig de vastgelegde gekalibreerde waarden zijn voor een toenemende omvang van het probleem. Dit is zowel relevant voor het kalibreren van het probleem in deze studie als het aanpassen van de parameters voor toepassing in praktijk wanneer het probleem groter of kleiner wordt afhankelijk van fluctuaties in de vraag en/of productiecapaciteit. De prestaties van de heuristiek met andere waarden voor populatiegrootte, selectiedruk en mutatiekans zijn quasi onveranderd. Enkel voor de populatie is een kleine verbetering van de doelfunctie merkbaar maar deze is niet significant. Omdat deze lichte verbetering niet gepaard gaat met significante verhoging van de rekentijd wordt toch besloten om deze parameterwaarde te verhogen naar 400. De impact van het aantal generaties (Figuur 7.10) en het aantal nakomelingen (Figuur 7.11) is echter veel sterker, maar dit gaat dan ook ten koste van een verhoogde rekentijd. Zowel de verbetering van de doelfunctie als de toename in rekentijd is veel sterker wanneer het aantal generaties wordt verhoogd in vergelijking met wanneer het aantal nakomelingen per generatie wordt verhoogd. 120%
10.000
9.000 100% 8.000 R² = 0,7973
7.000
6.000
60%
5.000
CPU Time
Objective function
80%
4.000 40% 3.000
2.000 20% 1.000
0%
0 50
100
150
200
250
# Generations Objective function
Computation time
Expon. (Objective function)
Figuur 7.10: Invloed van het aantal generaties op de rekentijd en de waarde van de doelfunctie bij toepassing op realistische probleemgrootte.
Een lineaire regressie wordt daarom gebruikt om aan de hand van de parameterwaarde de
Hoofdstuk 7. Genetische metaheuristiek
70
120%
10.000
9.000 100% R² = 0,715
8.000
7.000
6.000
60%
5.000
CPU Time
Objective function
80%
4.000 40% 3.000
2.000 20% 1.000
0%
0 200
300
400
500
# Offspring per generation Objective function
Computation time
Linear (Computation time)
Figuur 7.11: Invloed van het aantal nakomelingen op de rekentijd en de waarde van de doelfunctie bij toepassing op realistische probleemgrootte.
gemiddelde verbetering te gaan voorspellen voor wijzigingen van een van deze twee parameters. Aan de hand van deze resultaten kan dan worden bepaald wat het gemiddelde aantal seconden extra rekentijd zijn per procentuele verbetering van de doelfunctie ten opzichte van het referentiescenario. Door het verhogen van het aantal generaties kan gemiddeld genomen na 490 seconden een procent verbetering worden waargenomen in de doelfunctie. Echter bij het verhogen van het aantal nakomelingen per generatie zijn gemiddeld slechts 360 seconden nodig om een procent verbetering te brengen in de doelfunctie. Op basis hiervan wordt besloten om het aantal generaties op 50 te houden en het aantal nakomelingen per generatie te verhogen indien meer tijd of rekencapaciteit beschikbaar is.
7.2.7
Overzicht van parameters
Tabel 7.2 vat de waarden voor de parameters die in de bovenstaande secties werden bepaald nogmaals samen.
Hoofdstuk 7. Genetische metaheuristiek
Populatie omvang Pressure % Kans op mutatie
71
>400 10% - 20% 40%
Aantal generaties Aantal nakomelingen
50 ≥ 500
Tabel 7.2: Samenvatting van de waardes van de parameters zoals bepaald door de kalibratie.
7.3 7.3.1
Uitbreiding van de functionaliteiten Enkelvoudige heuristieken
Enkele oplossingen worden ook opgesteld aan de hand van zeer eenvoudige heuristieken. Deze heuristieken bestaan telkens uit twee delen en deze twee delen worden ook op verschillende manieren gecombineerd. Het eerste deel van de heuristieken verdeelt de jobs over de verschillende machines. De verschillende methodes die hiervoor gebruikt worden zijn hieronder in detail omschreven. Het tweede deel van de heuristieken gaat op zoek naar de beste volgorde van de jobs op de machine die ze toegewezen kregen. De verschillende methodes die hiervoor worden gebruikt worden eveneens hieronder beschreven. Selectie van machine Kortste productietijd Deze heuristiek overloopt voor elke job alle machines en plant de job in op de machine die de korste productietijd heeft. Indien er meerdere machines zijn met dezelfde productietijd wordt de job toegekend aan de machine waar op dat moment de laagste werklast aan toegewezen is. Gelijkmatige verdeling van werklast Hier wordt voor alle jobs sequentieel nagegaan welke machines de job kunnen produceren, uit deze groep machines wordt vervolgens de machine geselecteerd met de laagste werklast. Indien meerdere machines exact dezelfde werklast hebben wordt er willekeurig een machine geselecteerd uit deze machines. (Dit komt vooral voor aan het begin van het algoritme, wanneer op vele machines nog geen enkele job ingepland staat.)
Hoofdstuk 7. Genetische metaheuristiek
72
Dichtstbijzijnde klant Voor alle jobs wordt op basis van de vaste kosten geassocieerd met laattijdige levering bepaald welke machine het dichtst bij de klant gelegen is. Ook hier wordt eerst gecontroleerd of de machine in staat is om het order te produceren. Indien meerdere machines op dezelfde afstand van de klant staan wordt een willekeurige machine gekozen uit de machines die het dichtste bij de klant staan. Uitbalanceren van due dates
Deze heuristiek probeert dringende en minder dringende
jobs evenredig te verdelen over de verschillende machines. Dit gebeurt bij benadering door de som van de verschillende deadlines gelijkmatig te verdelen over de verschillende machines. Volgorde van jobs Vroegst beloofde datum Deze heuristiek gaat per machine de jobs rangschikken op de datum waarop de job beloofd werd aan de klant. Logischerwijs worden de jobs die het vroegst aan de klant beloofd werden het eerst geplaatst op de machine. Kortste productietijd Per machine worden de jobs gerangschikt op productietijd op de desbetreffende machine. De jobs met de kortste productietijd worden het eerst geplaatst om op die manier zo veel mogelijk orders binnen zo kort mogelijke tijd te kunnen afronden. Langste productietijd
Deze heuristiek doet het tegenovergestelde van de voorgaande, het
idee is hier om zeer laattijdige leveringen van orders te voorkomen door de langst durende orders als eerste te produceren. Kleinste verschil tussen release en beloofde datum
De berekeningen in deze heuristiek
zijn iets geavanceerder dan bij de bovenstaande. Voor alle jobs op de machine wordt het verschil gemaakt tussen de datum waarop de job beloofd wordt en de datum waarop de job wordt gereleased in het systeem (Dj − Rj ). De jobs waarvan het verschil tussen deze twee waardes het kleinste is worden als eerste geplaatst op de machines. Resultaten van verschillende combinaties De kwaliteit van deze verschillende single pass heuristieken wordt beoordeeld met behulp van de veertig parametersets die reeds gebruikt werden in sectie 5.2.1. De waarde van de
Hoofdstuk 7. Genetische metaheuristiek
73
doelfunctie wordt berekend voor elk schema dat kan worden gegenereerd met een combinatie van de bovenstaande heuristieken. Dit komt neer op 16 oplossingen per parameterset. Deze waarde wordt dan berekend met de optimale resultaten die bekomen werden door het vereenvoudigd optimalisatiemodel. Parametersets waarvoor geen optimale oplossing kon worden gevonden binnen redelijke tijd, en dus enkel bij benadering een oplossing hebben op basis van het vereenvoudigd optimalisatiemodel, werden uit deze vergelijkende studie weggelaten. Eerst worden de resultaten voor de verschillende machine selectie mechanismen besproken. Daarna worden de resultaten voor de volgorde heuristieken besproken. Ten slotte worden de beste combinaties van beide heuristieken besproken. De cijfers die gebruikt worden voor het vergelijken van de kwaliteit van de oplossingen van de heuristieken worden berekend door eerst de relatieve afwijking te berekenen in vergelijking tot de oplossing bekomen uit het vereenvoudigd optimalisatiemodel. Vervolgens worden deze relatieve afwijkingen gesommeerd over alle parametersets, om zo een representatief beeld weer te geven van de kwaliteit van de verschillende heuristieken over verschillende probleemgroottes. Average of Cumulative relative deviation
4500%
4000%
3500%
3000%
2500%
2000%
1500%
1000%
500%
0% M1
M2
M3
M4
Machine Sel
Figuur 7.12: Gemiddelde afwijking ten opzichte van de optimale oplossing voor de verschillende machine selectie mechanismen voor 40 willekeurige parametersets.
Op grafiek 7.12 wordt duidelijk dat de tweede en vierde methode voor het selecteren van
Hoofdstuk 7. Genetische metaheuristiek
74
machines duidelijk beter presteren dan de andere twee. Dit betekent dat het uitbalanceren van de belasting over de verschillende machines het beste presteert, op de voet gevolgd door de heuristiek die de due dates uitbalanceert over de verschillende machines heen. Average of Cumulative relative deviation
8000%
7000%
6000%
5000%
4000%
3000%
2000%
1000%
0% S1
S2
S3
S4
Sequence
Figuur 7.13: Gemiddelde afwijking ten opzichte van de optimale oplossing voor de verschillende volgorde toekennende mechanismen voor 40 willekeurige parametersets.
Een vergelijking van de prestaties van de verschillende volgorde mechanismen wordt gemaakt in grafiek 7.13. De vrij eenvoudige technieken die de vroegste beloofde datum en de kortste productietijd eerst plaatsen presteren hier het beste. De complexere vierde techniek presteert tegen de verwachtingen in slechter dan de andere heuristieken. De combinatie van beide soorten heuristieken in figuur 7.14 toont aan dat de prestaties van deze combinaties vooral worden bepaald door de heuristieken die de volgorde bepalen. Performantie van genetische heuristiek met opname van enkelvoudige heuristieken in startpopulatie De enkelvoudige heuristieken die hierboven werden beschreven worden nu gebruikt om een aantal oplossingen te genereren die opgenomen worden in de init¨ıele populatie van het genetisch algoritme. Deze heuristieken nemen iets meer rekentijd in beslag dan een willekeurige
Hoofdstuk 7. Genetische metaheuristiek
75
9000%
8000%
7000%
6000%
5000%
4000%
3000%
2000%
1000%
M4S4
M4S3
M4S2
M4S1
M3S4
M3S3
M3S2
M3S1
M2S4
M2S3
M2S2
M2S1
M1S4
M1S3
M1S2
M1S1
0%
Figuur 7.14: Gemiddelde cumulatieve afwijking ten opzichte van de optimale oplossing voor combinatie van machine en volgorde heuristieken voor 40 willekeurige parametersets.
generatie maar daar dit slechts eenmalig gebeurt bij het begin van het algoritme zorgt dit niet voor merkbare vertragingen.
7.3.2
Intelligente iteratie limieten voor het aantal generaties
Het bepalen van het juiste aantal generaties om de oplossing van het probleem zo optimaal mogelijk te maken is niet eenvoudig, vooral ook omdat het aantal nodige generaties sterk kan verschillen afhankelijk van de probleemgrootte. Een manier waarop dit kan opgelost worden is eenvoudigweg het gebruiken van informatie uit de voorgaande generaties om te beslissen of het nog nuttig is om een volgende generatielus uit te voeren. De twee meest bruikbare elementen om te bepalen of de voorige generaties nog effect hadden is het bepalen van de toename in de doelfunctie en het bepalen van het aantal nieuwe leden in de populatie. (Nieuwe leden worden enkel toegestaan in de populatie wanneer ze op zijn minst beter zijn dan het slechtste lid van de huidige populatie en geen identiek lid reeds voorkomt in de populatie, zie 7.1.6). Om te verzekeren dat het proces niet prematuur wordt stopgezet wordt hier gekozen voor de tweede methode, waarbij het aantal nieuwe leden dat wordt toegestaan in de populatie wordt
Hoofdstuk 7. Genetische metaheuristiek
76
bekeken. (Enkel de verbetering in het beste schema in beschouwing nemen zou kunnen leiden tot het vroegtijdig stopzetten van andere oplossingspaden die eveneens tot een beter resultaat kunnen leiden). Concreet wordt vooropgesteld dat wanneer gedurende de laatste g generaties geen nieuwe leden meer toegelaten werden tot de populatie de generatielus wordt stilgelegd, waarbij g een parameter is die vrij kan worden gekozen. Natuurlijk wordt er eveneens een bovengrens bepaald voor het aantal iteraties, om te voorkomen dat de berekeningen buitensporig veel tijd in beslag zouden nemen. Tests van deze methode wijzen erop dat vooral voor kleine problemen deze uitbreiding een merkbare verbetering geeft in rekentijd zonder de kwaliteit van de finale oplossing te verminderen.
7.4
Besluit
In dit hoofstuk werd het ontwerp, kalibratie en uitbreiding van een genetisch algoritme besproken. Verdere resultaten van deze meta-heuristiek zukken worden besproken in hoofdstuk 9.
Hoofdstuk 8
Simulated Annealing Metaheuristiek Een alternatief voor de genetische metaheuristiek die vertrekt van een populatie om een optimale oplossing te zoeken is de op lokaal zoeken gebaseerde simulated annealing metaheuristiek. Deze techniek is gebaseerd op het gedrag van atomen tijdens het afkoelingsproces van metalen. Hoe hoger de temperatuur van het metaal is, hoe vrijer de deeltjes in het metaal kunnen bewegen. Deze fysische reactie wordt gebruikt als analogie voor een metaheuristiek die afhankelijk van de temperatuur (een parameter van de heuristiek) een verschillende kans berekent voor het aanvaarden van minder goede oplossingen. De precieze werking van deze techniek wordt uitgebreid besproken in 8.1. De heuristiek die hier werd ontwikkeld is in grote lijnen gebaseerd op het werk van Celano et al. (2008) die voor een analoog probleem werkt met een simulated annealing heuristiek in combinatie met een simulatie ter evaluatie van de kwaliteit van de schema’s. Daarnaast werden zowel willekeurige lokale zoekbewegingen als intelligente zoekbewegingen gebruikt, gebaseerd op het werk van Piersma & van Dijk (1996).
8.1
Werking van simulated annealing metaheuristiek
Deze sectie geeft een gedetailleerd overzicht van de werking van deze metaheuristiek. Eerst wordt gestart met een globaal overzicht van de verschillende stappen van de heuristiek en vervolgens wordt op relevante stappen dieper ingegaan.
77
Hoofdstuk 8. Simulated Annealing Metaheuristiek
8.1.1
78
Basisprincipes van simulated annealing
De simulated annealing metaheuristiek start van een enkele oplossing. Deze oplossing kan willekeurig worden gegenereerd of gebaseerd zijn op enkelvoudige heuristieken. In vergelijking met het genetisch algoritme is de zoekruimte dus minder divers maar dit wordt gecompenseerd door het aanvaarden - onder voorwaarden - van oplossingen die een minder goed resultaat hebben dan de huidige oplossing. Dit betekent dat de enkele oplossing steeds verder evolueert tot na een bepaald aantal iteraties het algoritme wordt stopgezet. Tijdens elke iteratie van de simulated annealing heuristiek wordt de huidige beste oplossing vergeleken met een alternatieve oplossing die bekomen is door middel van lokaal zoeken. Indien deze nieuwe oplossing beter is dan de huidige beste oplossing wordt deze vervangen door de nieuwe betere oplossing.
P (accept) = exp(−
CU RREN T − BEST ) Tc
(8.1)
Wanneer de nieuwe oplossing slechter is dan de huidige beste oplossing is er een kans dat deze alsnog geaccepteerd wordt, afhankelijk van een kans die wordt berekend op basis van de huidige temperatuur en het verschil tussen beide oplossingen. Vergelijking 8.1 berekent de kans dat een slechter schema toch wordt geaccepteerd op basis van deze parameters. Dit betekent dat wanneer de temperatuur hoger is een slechtere oplossing een hogere kans heeft om toch geaccepteerd te worden. Naarmate het algoritme loopt wordt de temperatuur steeds lager, opdat de oplossing optimaal zou convergeren. Naast de temperatuur zorgt een kleiner verschil (= een kleiner negatief effect op de doelfunctie) ook voor een grotere kans op accepteren van de slechtere oplossing. De belangrijkste parameters van de simulated annealing heuristiek worden samengevat in tabel 8.1. De procedure gaat als volgt in zijn werk: Stap 1 De waarden van de verschillende parameters worden ingelezen en het op te lossen probleem wordt ingeladen. Stap 2 Een startoplossing (s) wordt willekeurig of op basis van de enkelvoudige heuristieken uit 7.3.1 gegenereerd. De waarde van de doelfunctie voor deze oplossing wordt ge¨ınitializeerd als huidige beste oplossing. (BEST)
Hoofdstuk 8. Simulated Annealing Metaheuristiek
Parameter
79
Beschrijving
Ts
Begintemperatuur waarvan de afkoeling gestart wordt.
Tf
Eindtemperatuur waarbij de heuristiek wordt stopgezet.
α
Ratio die snelheid van afkoelen per generatie weergeeft.
Nmax
Maximaal aantal lokale bewegingen voor temperatuurverlaging.
Prand
Kans op willekeurige lokale zoekbewegingen.
Pswap
Kans op swap operatie in local search. Tabel 8.1: Parameters van simulated annealing heuristiek
Stap 3 De huidige temperatuur wordt gelijkgesteld aan de starttemperatuur (Tc = Ts .) Stap 4 De teller voor lokale bewegingen wordt ge¨ınitialiseerd (N = 0). Stap 5 Genereer een willekeurig getal y tussen 0 en 1; als y < Prand gebruik dan een willekeurige methode voor het lokaal zoeken, gebruik anders de inteligente methode voor lokaal zoeken. Stap 6 Genereer een willekeurig getal z tussen 0 en 1; als z < Pswap genereer dan een nieuw schema (s0 ) met behulp van een swap operatie, anders wordt er een nieuw schema gegenereerd door verplaatsing. De doelfunctie van dit nieuwe schema wordt berekend en opgeslaan als CURRENT; N = N+1. Stap 7 Als CU RREN T < BEST stel dan s = s0 en stel N = Nmax en ga naar stap 10 Stap 8 Genereer een willekeurig getal x tussen 0 en 1. Stap 9 Als x ≤ exp(− CU RRENTcT −BEST ) stel dan s = s’ en stel N = Nmax . Stap 10 Als N < Nmax ga dan naar stap 5, anders stel Tc = αTc . Stap 11 Als Tc > Tf ga naar stap 4. Stap 12 STOP, informatie wordt opgeslaan.
8.1.2
Methodes voor lokaal zoeken binnen de simulated annealing heuristiek
Voor het lokaal zoeken worden twee verschillende methodes gebruikt: verwisselen en verplaasten. Verplaatsen houdt in dat een enkele job een andere plaats krijgt toegewezen binnen de
Hoofdstuk 8. Simulated Annealing Metaheuristiek
80
planning. Verwisselen houdt in dat twee jobs met elkaar van plaats gaan wisselen. Welke jobs hiervoor worden geselecteerd kan op verschillende manieren worden bepaald. Deze selectiemethoden kunnen zowel volledig willekeurig als op een bepaalde manier berekend zijn. De intelligente manier om jobs te kiezen werd hier gebaseerd op de methode van Piersma & van Dijk (1996). Deze methode gaat op zoek naar de zwaarst belaste machines en probeert deze bij voorkeur te ontlasten. (Onder het zwaar belast zijn van een machine wordt hier verstaan een machine waarbij de gecumuleerde bijdrage van de aanwezige jobs tot de doelfunctie relatief groot is. Doordat de doelfunctie in Piersma & van Dijk (1996) anders was dan de doelfunctie die hier gebruikt wordt is de implementatie van dit principe anders.) Willekeurig lokaal zoeken Deze sectie geeft wat meer duiding bij de werking van het willekeurig lokaal zoeken. Veronderstel de situatie zoals in figuur 8.1 waarbij 15 jobs op 5 machines ingepland dienen te worden. Het willekeurig lokaal zoeken gebeurt door lijsten op te stellen voor alle mogelijke verplaatsingen en verwisselingen van jobs. Deze lijsten worden vervolgens willekeurig dooreengeschud. Tijdens het uitvoeren van de iteraties wordt dan afhankelijk van de waarde van Pswap de volgende beweging geselecteerd uit de lijst met wissels ofwel de lijst met verplaatsingen. Als we kijken naar de verplaatsingen die kunnen gemaakt worden in de situatie van figuur 8.1 zien we dat elke job naar 19 verschillende locaties kan verplaatst worden. Dit betekent dat in het totaal 285 verplaatsingen mogelijk zijn (= 19 verplaatsingsmogelijkheden * 15 jobs).
Hoofdstuk 8. Simulated Annealing Metaheuristiek
81
Figuur 8.1: Voorbeeldschema ter illustratie van simulated annealing metaheuristiek
Voor het mogelijke aantal swaps (= twee jobs die met elkaar van plaats verwisselen) wordt een combinatie van 2 uit een populatie van 15 jobs berekend om alle mogelijke wissels te bekomen. Dit komt neer op 105 (=
15! 2!13! )
mogelijke onderlinge wissels.
Intelligent lokaal zoeken Voor het intelligent zoeken wordt uitgegaan van dezelfde mogelijke verplaatsingen en swaps, maar wordt berekend welke bewegingen waarschijnlijk het beste effect op de doelfunctie zullen hebben. Dit kan leiden tot snellere verbeteringen van de kwaliteit van het schema. Een nadeel van deze methode is natuurlijk dat hiervoor meer rekentijd nodig is. Daarnaast is er ook nog een belangrijker probleem, namelijk dat deze methode gevoeliger is om te blijven steken in lokale optima. De oorzaak van dit probleem is tweeledig. Allereerst is het minder waarschijnlijk dat slechtere bewegingen worden geselecteerd, wat betekent dat minder opwaartse bewegingen zullen worden uitgevoerd, en lokale optima dus automatisch worden bevoordeeld. Ten tweede bestaat de mogelijkheid bij zeer grote problemen dat bepaalde bewegingen nooit geprobeerd zullen worden omdat deze verwacht worden een negatief effect te hebben op de doelfunctie. Bij een probleem van realistische groote dat 300 jobs op 32 machines moet gaan plannen zijn er in totaal 144.450 wissels en verplaatsingen mogelijk, het is duidelijk dat wanneer enkel inteligent zoeken zal worden gebruikt een groot deel van deze bewegingen geen
Hoofdstuk 8. Simulated Annealing Metaheuristiek
82
enkele kans heeft om ooit aan bod te komen doordat het aantal uitgeprobeerde bewegingen wordt beperkt door Nmax . Dit is de reden waarom hier werd geopteerd voor een aanpak die beide methodes combineert. Hieronder wordt besproken hoe het sorteren van deze twee verzamelingen precies gebeurt. (Zoals eerder vermeld zijn deze methodes gebaseerd op Piersma & van Dijk (1996).) De paragrafen hieronder bespreken hoe de lokale zoekbewegingen (wissels en verplaatsingen) worden gesorteerd voor gebruik in intelligente zoekmethodes.
Sorteren van mogelijke verplaatsingen De eerste stap in het sorteren van de mogelijke verplaatsingen van de jobs is het berekenen van de contributie tot de doelfunctie van de verschillende machines. De redenering is hier dat bewegingen van jobs op deze machines de grootste kans hebben om de doelfunctie in het algemeen te gaan verbeteren. De machine met de grootste contributie tot de doelfunctie wordt geselecteerd, en voor alle jobs op deze machine worden de effici¨enties berekend van de mogelijke verplaatsingen naar andere machines. Deze worden berekend aan de hand van formule 8.2. In deze formule is mj de minimale productietijd van job j (= de productietijd op de meest efficiente machine die in de dataset is opgenomen.) en Pjm de productietijd van job j op machine m. Op deze manier zal de effici¨entie dus steeds een getal tussen 0 en 1 zijn. (Noot: Bij het intelligent zoeken worden dus enkel de bewegingen van de zwaarst belaste machine naar andere machines in beschouwing genomen.) Voor het voorbeeld in figuur 8.1 wordt deze berekening ge¨ıllustreerd met behulp van tabel 8.2 die de productietijden bevat (een waarde M betekent in deze tabel dat productie van deze job niet mogelijk is op deze machine) en tabel 8.3 waarin de berekende effici¨enties opgelijst staan.
ef f (j, m) =
mj Pjm
(8.2)
In dit voorbeeld heeft de eerste machine de grootste contributie tot de doelfunctie. De verplaatstingen van jobs die op deze machine ingeplaatst zijn worden dus in beschouwing genomen. De effici¨enties voor de verschillende mogelijke verplaatsingen zijn in het vet aangeduid in tabel 8.2.
Hoofdstuk 8. Simulated Annealing Metaheuristiek
83
Machine number Job
1
2
3
4
5
1
0,70
0,60
M
0,60
M
2
6,10
5,24
M
5,24
6,10
3
10,20
8,76
10,20
8,76
M
4
7,46
6,40
M
6,40
M
5
10,72
11,48
M
10,72
11,48
6
12,92
11,10
M
11,10
M
7
6,12
5,26
M
5,26
M
8
5,12
4,39
4,39
4,39
M
9
3,41
2,93
M
2,93
M
10
8,03
6,90
6,04
6,90
6,44
11
0,54
0,58
M
0,54
M
12
15,75
16,87
M
15,75
15,75
13
6,10
5,24
M
5,24
M
14
15,13
12,13
11,38
12,99
M
15
0,72
0,62
0,72
0,62
M
Tabel 8.2: Productietijden op de verschillende machines.
Hoofdstuk 8. Simulated Annealing Metaheuristiek
84
Machine number Job
1
2
3
4
5
1
0,86
1,00
0,00
1,00
0,00
2
0,86
1,00
0,01
1,00
0,86
3
0,86
1,00
0,86
1,00
0,01
4
0,86
1,00
0,01
1,00
0,01
5
1,00
0,93
0,01
1,00
0,93
6
0,86
1,00
0,01
1,00
0,01
7
0,86
1,00
0,01
1,00
0,01
8
0,86
1,00
1,00
1,00
0,00
9
0,86
1,00
0,00
1,00
0,00
10
0,75
0,88
1,00
0,88
0,94
11
1,00
0,93
0,00
1,00
0,00
12
1,00
0,93
0,02
1,00
1,00
13
0,86
1,00
0,01
1,00
0,01
14
0,75
0,94
1,00
0,88
0,01
15
0,86
1,00
0,86
1,00
0,00
Tabel 8.3: Effici¨enties op de verschillende machines.
Hoofdstuk 8. Simulated Annealing Metaheuristiek
85
Let wel dat er veel meer verplaatsingen mogelijk zijn dan er effici¨entiescores berekend worden. Alle mogelijke verplaatsingen worden gesorteerd en de verplaatsingen met de hoogste effici¨entiescore zullen dan als eerste in aanmerking komen voor uitvoering. Wanneer bijvoorbeeld de verplaatsing van job 11 naar de eerste positie van machine 2 wordt vergeleken met de verplaatsing van job 13 naar machine 5 zien we dat de verplaatsing van job 11 eerst zal worden geprobeerd, de verplaatsing van job 13. Deze rangschikkingen zijn zowel afhankelijk van de job die verplaatst wordt als de machine naar waar verplaatst wordt. Deze scores zijn echter niet afhankelijk van de plaats waarnaar de job verplaatst wordt. De reden hiervoor is dat het effect van de verplaatsing niet eenduidig is omdat jobs heel uiteenlopende kenmerken hebben (deadline, productietijd...) waardoor het effect op de algemene doelfunctie niet altijd voor de hand liggend is.
Sorteren van mogelijke swaps
Voor het sorteren van de swaps wordt eveneens gebruik
gemaakt van de effici¨enties om de mogelijkheden te gaan sorteren. Hier wordt echter niet alleen de machine met de hoogste contributie tot de doelfunctie in beschouwing genomen, maar worden alle mogelijke wissels overlopen. Hieronder worden de verschillende stappen voor het sorteren van de swaps uitgeschreven. Stap 1 Sorteer alle machines op afnemende volgorde van contributie aan de doelfunctie: Machines = { m1 ,m2 ,..., mN M } waarbij m1 de hoogste en mN M de laagste contributie heeft. Stap 2 Stel dan mv = m1 en mn = mN M . stap 3 Sorteer de jobs op Mv in afnemende volgorde van effici¨entie op Mn en sorteer de jobs op Mn in afnemende volgorde van effici¨entie op Mv . De geordende verzamelingen van jobs krijgen als naam respectievelijk: Jv = {jv1 , jv2 , ..., jvN V } en Jn = {jn1 , jn2 , ..., jnN N }. Stap 4 Voeg alle mogelijke wissels tussen deze twee machines toe aan de geordende verzameling van wissels S. (In de eerste plaats wordt de volgorde op machine mv gerespecteerd, daarna de volgorde op machine mn . Stap 5 Stel Mn = mn−1 , als Mn 6= Mv ga naar stap 3, anders ga naar stap 6. Stap 6 Als Mv ≥ mN M −1 stel dan Mv = Mv+1 en mn = mN M en ga naar stap 3, anders ga naar stap 7.
Hoofdstuk 8. Simulated Annealing Metaheuristiek
86
Stap 7 STOP
De stappen die hierboven worden ge¨ıllustreerd in figuur 8.2. Het eerste luik van de figuur toont de beginsituatie analoog aan die uit figuur 8.1. In de eerste stap worden de machines gerangschikt in volgorde van contributie aan de doelfunctie, van hoog naar laag. Dan worden in de tweede stap de eerste en de laatste machine geselecteerd om de mogelijke wissels toe te voegen aan een geordende verzameling van wissels. De volgorde van de jobs op deze machines wordt tevens aangepast om de effici¨entie van productie op de andere machine te weerspiegelen. Ter illustratie worden in tabel 8.3 de toegevoegde wissels weergegeven (in de volgorde waarin deze werden toegevoegd). Vervolgens worden de wissels met machines met minder hoge contributies tot de doelfunctie overlopen in stap 5, waarna in stap zes de machine met de op een na hoogste contributie tot de doelfunctie de plaats van de machine met de hoogste contributie zal innemen. Dit proces wordt verdergezet tot de wissels tussen de wissels tussen de machines met de op een na laagste en de laagste contributie toegevoegd zijn aan de geordende lijst met wissels.
Hoofdstuk 8. Simulated Annealing Metaheuristiek
Stap 0 Machine 1 Machine 2 Machine 3 Machine 4 Machine 5
1 2 3 4 5
6 7 8 9 10
11 12 14
87
13 15
Stap 1 Machine 1 Machine 4 Machine 5 Machine 3 Machine 2
1 4 5 3 2
6 9 10 8 7
11
13
14 12
15
Stap 3 Machine 1 Machine 4 Machine 5 Machine 3 Machine 2
1 4 5 3 12
6 9 10 8 2
11
13
14 7
15
Stap 4 & 5 Machine 1 Machine 4 Machine 5 Machine 3 Machine 2
1 4 5 3 12
6 9 10 8 2
11
13
15 7
14
Stap 6 Machine 1 Machine 4 Machine 5 Machine 3 Machine 2
1 4 5 3 12
6 9 10 8 2
11
13
15 7
14
Doelfunctie 3.922 18 1.441 1.705 1.476 Doelfunctie 3.922 1.705 1.476 1.441 18 Doelfunctie 3.922 1.705 1.476 1.441 18 Doelfunctie 3.922 1.705 1.476 1.441 18 Doelfunctie 3.922 1.705 1.476 1.441 18
Figuur 8.2: Bepalen van de volgorde van de wissels bij intelligent zoeken.
8.2
Kalibratie van parameters voor simulated annealing metaheuristiek
In deze sectie wordt overgegaan tot de kalibratie van de parameters van de simulated annealing metaheuristiek. Hiervoor wordt eerst een manuele inschatting gedaan van de grootteorde van de verschillende parameters in sectie 8.2.1. Vervolgens worden de parameters een voor een gekalibreerd.
Hoofdstuk 8. Simulated Annealing Metaheuristiek
88
Nr
m1
m2
1
1
12
2
1
2
3
1
7
4
6
12
5
6
2
6
6
7
7
11
12
8
11
2
9
11
7
10
13
12
11
13
2
12
13
7
Tabel 8.4: wissels toegevoegd aan geordende lijst in stap 4 uit figuur 8.2
De grafische voorstelling van de kalibratie gebeurt steeds aan de hand van figuren zoals figuur 8.4. Deze voorstelling is analoog aan de voorstelling in sectie 7.2 waar de parameters voor de genetische metaheuristiek werden gekalibreerd.
8.2.1
Schatting van optimale parameterwaarden
Een eerste schatting van de parameterwaarden kan worden gemaakt door het aantal iteraties en de gemiddelde kansen voor het accepteren van een slechte oplossing in deze verschillende iteraties te gaan berekenen. De basis voor het schatten van deze kansen is de gemiddelde impact van een lokale zoekbeweging op de waarde van de doelfunctie. Op basis van enkele runs werd deze waarde vastegelegd op 15.000. Het exacte aantal iteraties dat nodig zal zijn kan worden berekend met de formule: logα ( TTfs ). De geschatte startwaarden op basis van de bovenstaande assumpties worden samengevat in tabel 8.5. Deze keuze werd gemaakt omdat deze waarden een start met hoge kansen op accepteren van slechtere oplossingen garanderen, evenals een langere periode waarin de kans op het accepteren van een slechtere oplossing quasi nihil is. Ook is het aantal iteraties beperkt, wat een kortere rekentijd oplevert voor de eerste tests.
Hoofdstuk 8. Simulated Annealing Metaheuristiek
Ts
50.000
Tf
100
α
95%
Nmax
500
Prand
50%
Pswap
50%
89
Tabel 8.5: Startwaarden voor kalibratie van simulated annealing heuristiek. 80%
60.000
70% 50.000
40.000
Temperature
50%
40%
30.000
30% 20.000
Probability of acceptance
60%
20% 10.000 10%
0
0% 1
5
9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89 93 97 101105109113117121
Temperature
P(accept)
Figuur 8.3: Gemiddelde kans op accepteren van nieuwe oplossing met eerste schattingen van parameters.
Deze gemiddelde verwachtingen voor deze waarden worden grafisch voorgesteld in figuur 8.3.
8.2.2
Starttemperatuur
De waarde van de starttemperatuur wordt gevarieerd van 10.000 tot 250.000, de resultaten worden samengevat in grafiek 8.4. De resultaten van de kalibratie voor deze parameter tonen aan dat zolang deze parameter binnen de juiste grootteorde wordt gehouden de impact op het resultaat vrij gering is. Ook met betrekking tot rekentijd blijft de impact hier beperkt. Natuurlijk zal een hogere starttemperatuur tot gevolg hebben dat een groter aantal iteraties wordt overlopen maar hier blijkt vooral de waarde voor α belangrijk. Het veiligste lijkt hier om de waarde eerder aan de bovenzijde van het spectrum vast te leggen.
Hoofdstuk 8. Simulated Annealing Metaheuristiek
90
Objective CPU Time
500
3.000.000
450 400 350
2.000.000
300 250
1.500.000
200 1.000.000
CPU time (sec)
Objective function value
2.500.000
150 100
500.000
50 0
0 10.000
20.000
30.000
40.000
50.000
100.000
150.000
200.000
250.000
Starting temperature Values
Objective
CPU Time
T_start
Figuur 8.4: Invloed van de starttemperatuur op de doelfunctie en de rekentijd.
Dit omdat alle waarden onder de startwaarde sowieso worden overlopen in de iteraties tot de eindtemperatuur. De waarde wordt daarom vastgelegd op 200.000.
8.2.3
Eindtemperatuur
Voor de eindtemperatuur worden de waarden gevarieerd van 101 tot 1. De resultaten worden weergegeven in grafiek 8.5. Analoog met de starttemperatuur kan hier worden geconcludeerd dat de impact van deze parameter op het resultaat gering is. De impact op de rekentijd wordt ook hier vooral gedomineerd door de parameter α. Wel is zichtbaar dat voor een waarde van 1 de kwaliteit van de oplossing een merkbare verbetering toont, maar dit is te wijten aan het grotere aantal iteraties. (Het grotere aantal iteraties zorgt er ook metteen voor dat de rekentijd gevoelig groter wordt) Omdat dit aantal iteraties eenvoudiger kan worden aangepast met de parameter α wordt besloten om de waarde van deze parameter aan de bovenzijde van het spectrum te houden opdat deze de rekentijd zo min mogelijk zou be¨ınvloeden. Deze parameter wordt op 100 vastgelegd, zolang hier binnen de juiste grootteorde wordt gebleven levert deze parameter geen problemen op.
Hoofdstuk 8. Simulated Annealing Metaheuristiek
91
Objective CPU Time
500
3.500.000
450 3.000.000 400 350 300
2.000.000
250 1.500.000
200
Computation time (sec)
Objective function value
2.500.000
150
1.000.000
100 500.000 50 -
1
11
21
31
41
51
61
71
81
91
101
Ending temperature Values
Objective
CPU Time
T_stop
Figuur 8.5: Invloed van de eindtemperatuur op de doelfunctie en de rekentijd.
8.2.4
Alfa
Uit figuur 8.6 blijkt dat de invloed van de alfa parameter op de kwaliteit van de oplossing en de rekentijd is veel sterker is. De kwaliteit van oplossing neemt quasi evenredig toe met toenemende waarde van deze parameter. Natuurlijk neemt ook de oplossingstijd sterk toe waneer deze parameter stijgt omdat onmiddelijk een veel groter aantal iteraties nodig zijn alvorens de eindtemperatuur wordt bereikt. Het lijkt hier aangewezen om de waarde van deze parameter zo hoog mogelijk te laten oplopen zonder de rekentijd buitensporige proporties te laten aannemen. Hoe groot een redelijke rekentijd is hangt natuurlijk af van de snelheid van de computers en de beschikbare tijd om het model uit te gaan rekenen. Voorlopig wordt de waarde van deze parameter ingesteld op 99.5% omdat dit het punt lijkt te zijn waar het grootste deel van het opwaarts potentieel voor de doelfunctie al is bereikt, en de oplossingstijd nog niet exponentieel begint toe te nemen. Om de prestaties van de heuristiek te verbeteren kan deze waarde nog worden verhoogd indien gewenst.
Hoofdstuk 8. Simulated Annealing Metaheuristiek
92
Objective functionCPU Time
9.000
3.000.000
8.000 2.500.000
6.000
2.000.000
5.000 1.500.000 4.000
3.000
1.000.000
Computation time (sec)
Objective function value
7.000
2.000 500.000 1.000
0
0 95,0%
96,0%
97,0%
98,0%
99,0%
99,1%
99,2%
99,3%
99,4%
99,5%
99,6%
99,7%
99,8%
99,9%
Alpha parameter value Values
Objective function
CPU Time
Alpha
Figuur 8.6: Invloed van alfa op de doelfunctie en de rekentijd.
8.2.5
Maximaal aantal lokale zoekbewegingen
De resultaten van de kalibratie voor Mmax worden weergegeven in figuur 8.7. Hier blijkt dat deze parameter geen al te sterke impact heeft op zowel de resultaten als de rekentijd. Dit is contra-intu¨ıtief omdat de verwachting is dat de rekentijd zal toenemen wanneer het algoritme wordt toegelaten om meer zoekbewegingen uit te voeren per iteratie. Omdat deze impact op de rekentijd niet erg groot is kan deze waarde vrij hoog worden vastgelegd om eventuele restricties op de zoekruimte te gaan vermijden. Daarom wordt het maximale aantal bewegingen vastgelegd op 20.000.
8.2.6
Kans op willekeurig zoeken
Figuur 8.8 geeft de resultaten van de kalibratie van de kans op willekeurig zoeken weer. Afhankelijk van de waarde van deze parameter zal gemiddeld genomen een bepaalde fractie van de iteraties gebruik maken van het willekeurig zoekmechanisme of het intelligent zoekmechanisme. Een waarde van 25% voor deze parameter betekent dat gemiddeld gezien 25% van de iteraties gebruik zullen maken van willekeurig lokaal zoeken en in 75% van de iteraties van intelligent lokaal zoeken.
Hoofdstuk 8. Simulated Annealing Metaheuristiek
93
Average of Final objectiveAverage of CPU Time
1.800
800.000
1.600
700.000
1.400
1.200 500.000 1.000 400.000 800 300.000 600 200.000
Computation time (sec)
Objective function value
600.000
400
100.000
200
0
0 500
1500
2500
3500
4500
5000
10000
15000
20000
25000
30000
35000
40000
45000
50000
55000
60000
Maximum local search moves per iteration Values Average of Final objective
Average of CPU Time
M_iterations
Figuur 8.7: Invloed van van Mmax op doelfunctie en rekentijd.
Uit de resultaten blijkt duidelijk dat een hogere waarde voor deze variabele zowel betere resultaten als een lagere rekentijd geeft. De verlaging in rekentijd is logisch want voor het sorteren van de lokale zoekbewegingen is natuurlijk enige rekentijd nodig. Minder verwacht is dat de willekeurige zoekbewegingen gemiddeld blijkbaar beter presteren dan de intelligente zoekbewegingen. Het is duidelijk dat deze parameter best een waarde aan de bovenkant van het spectrum toegewezen krijgt. De waarde voor deze parameter wordt vastgelegd op 90%.
8.2.7
Kans op wissel
Figuur 8.9 geeft de resultaten weer van de kalibratie van de kans op wissel. Deze parameter bepaalt in welke fractie van de zoekbewegingen gebruik wordt gemaakt van een wissel of een verplaatsing. Een waarde van 90% voor deze parameter zou betekenen dat 90% van de zoekbewegingen een wissel zijn en 10% van de gevallen een verplaatsing van een job. De resultaten tonen duidelijk dat wanneer een hoger aandeel wissels wordt gebruikt in verhouding tot verplaatsingen die doelfunctie een hogere - en dus slechtere - waarde heeft. Dit kan verklaard worden doordat bij een wissel steeds een job zal terugkeren naar een zwaarder belaste machine en zo natuurlijk minder zal ontlast worden dan wanneer een job gewoon
Hoofdstuk 8. Simulated Annealing Metaheuristiek
94
3.500.000
1.400
3.000.000
1.200
2.500.000
1.000
2.000.000
800
1.500.000
600
1.000.000
400
500.000
200
Computation time (sec)
Objective function value
Objective CPU time
0
0 10,00%
20,00%
30,00%
40,00%
50,00%
60,00%
70,00%
80,00%
90,00%
Chance of random local search Values
Objective
CPU time
pRandom
Figuur 8.8: Invloed van kans op willekeurig zoeken op doelfunctie en rekentijd.
verdwijnt van de machine. Wel is er een toename - ongeveer 30% - in rekentijd waar te nemen wanneer het percentage wissels zeer laag wordt (10%). Dit kan worden verklaard doordat het algoritme meer in dezelfde groep van mogelijke verplaatsingen zal zoeken, waardoor op het einde van het algoritme langer zal moeten worden gezocht voordat een positieve beweging wordt gevonden. Het lijkt dus zinnig om deze parameter op 20% vast te stellen, op deze manier wordt een oplossing van betere kwaliteit bekomen en wordt terzelfdertijd de rekentijd ietwat ingeperkt.
8.2.8
Overzicht van gekalibreerde parameterwaarden voor simulated annealing
Tabel 8.6 toont een overzicht van de parameterwaarden die in de in de vorige secties werden bepaald.
8.3
Besluit
In dit hoofstuk werd een genetisch algoritme ontworpen en gekalibreerd voor optimale prestaties bij het oplossen van het bedrijfsspecifieke probleem. Resultaten en performatie van dit
Hoofdstuk 8. Simulated Annealing Metaheuristiek
95
Objective CPU Time 3.500.000
1.000
900 3.000.000 800
700
600 2.000.000 500 1.500.000 400
Computation time (sec)
Objective function value
2.500.000
300
1.000.000
200 500.000 100
0
0 10,00%
20,00%
30,00%
40,00%
50,00%
60,00%
70,00%
80,00%
90,00%
Probability of using swap move Values Objective
CPU Time
pSwap
Figuur 8.9: Invloed van kans op wissel op doelfunctie en rekentijd.
Ts
200.000
Tf
100
α
99.5%
Nmax
20.000
Prand
90%
Pswap
20%
Tabel 8.6: Overzicht van de gekalibreerde waarden voor de simulated annealing heuristiek.
algoritme worden verder besproken in hoofdstuk 9.
Deel III
Resultaten en uitbreidingen
96
Hoofdstuk 9
Performantie van heuristieken In dit hoofdstuk worden de heuristieken die in hoofstukken 7 en 8 zijn ontwikkeld en gekalibreerd uitgebreid getest. Het testen gebeurt op twee manieren, enerzijds worden de beide heuristieken vergeleken met de exacte oplossingsmethode en anderzijds worden deze getest op problemen van realistische groote (zes weken productie op 32 machines).
9.1
Vergelijking met exacte oplossingstechnieken
Grafiek 9.1 toont de gemiddelde relatieve performantie van de twee heuristieken en de exacte oplossingsmethode wanneer deze op eenzelfde probleem worden uitgetest (deze resultaten zijn gebaseerd op 25 willekeurig gegenereerde parametersets). De eerste conclusie is dat alle drie de technieken zeer grote verbeteringen kunnen aanbrengen wanneer deze worden vergeleken met de huidige planningstechnieken. Deze laatste worden benaderd met behulp van enkelvoudige heuristieken. De maximale verbetering van 70% wordt hier aangegeven door de resultaten van de exacte optimalisatie. Natuurlijk heeft deze methode de langste oplossingstijd nodig. Door het zeer grote verschil in oplossingstijden worden de oplossingstijden hier op een logaritmische schaal weergegeven, andere assen zijn niet op logaritmische schaal. Het genetisch algoritme presteert merkbaar minder goed dan het simulated annealing algoritme wanneer beide heuristieken worden vergeleken. Deze prestatie is zowel op vlak van eindresultaat als op vlak van rekentijd inferieur aan die van de simulated annealing heuristiek. Het simulated annealing algoritme brengt 36% meer verbetering wanneer het wordt vergeleken met de oplossing van de genetische heuristiek. 97
Hoofdstuk 9. Performantie van heuristieken
98
10.000
70.000
60.000
Objective function value
40.000 100 30.000
20.000
CPU time (seconds, LOG SCALE)
1.000
50.000
10
10.000
1
0
Single pass heuristic
Exact optimization Genetic Algorithm Solution method
Objective function
Simulated Annealing
Computation time
Figuur 9.1: Performantie van exacte en meta heuristische oplossingsmethodes getest op probleem dat 15 jobs inplant op 4 machines.
Ook voor de rekentijd is het simulated annealing algoritme superieur, de rekentijd van het genetisch algortime is hier een factor 12 hoger dan de rekentijd voor het simulated annealing algoritme. Doordat deze resultaten gebaseerd zijn op problemen van beperkte grootte is het niet mogelijk om deze gegevens zomaar te aanvaarden. Daarom wordt in de volgende sectie meer aandacht besteed aan de toepassing van beide algoritmen op problemen die eenzelfde omvang hebben als de problemen die in praktijk opgelost moeten worden. Wel geeft deze informatie aan dat zelfs voor zeer beperkte schaal de complexiteit van deze problemen niet onderschat mag worden.
9.2 9.2.1
Toepassing op realistische probleemgrootte Performantie met standaardparameters uit kalibratie
Voor het testen van de prestaties van de heuristieken op problemen van realistische grootte werden 20 verschillende parametersets gegenereerd die elk 300 jobs bevatten die ingepland moeten worden op 32 machines verdeeld over twee productielocaties. Deze problemen zijn
Hoofdstuk 9. Performantie van heuristieken
99
veel te complex om met exacte methodes te worden opgelost. De best presterende enkelvoudige heuristiek wordt hier beschouwd als een goede oplossing volgens de huidige ad hoc planningsmethode. In realiteit zal de planning die nu wordt gebruikt vaak zelfs deze kwaliteit niet bereiken. De genetische meta-heuristiek toont gemiddeld over alle oplossingen een verbetering van 43% in de doelfunctie ten opzichte van de enkelvoudige heuristieken. Hiervoor heeft deze methode (met de standaardparameters zoals die werden gekalibreerd in sectie 7.2.7) gemiddeld 3.748 seconden of dus ongeveer 1 uur nodig. Simulated Annealing vertoont echter een nog drastischere verbetering van 84% ten opzichte van de enkelvoudige heuristieken. Daarenboven is deze heuristiek (met de gekalibreerde waarden uit sectie 8.2.8) sneller en heeft deze gemiddeld slechts 1.206 seconden of dus ongeveer 20 minuten nodig. Hieruit blijkt dat het simulated annealing algoritme de duidelijke winnaar is met oplossingen die gemiddeld gezien 72% beter zijn van kwaliteit dan de oplossingen uit het genetisch algoritme. Daarnaast kan ook worden gesteld dat het verschil in performantie verder toeneemt wanneer de omvang van het probleem toeneemt.
9.2.2
Performantie met hogere toegelaten rekentijd
Bij de tweede test werd opnieuw een parameterset van realistische grootte gegenereerd. Opnieuw worden beide algoritmen getest op deze parameterset, maar deze keer wordt bij elk algoritme een variabele enkele malen van waarde veranderd zodat het algoritme een langere rekentijd heeft om het probleem op te lossen. Voor het genetisch algoritme is dit het aantal nakomelingen. Zoals eerder besproken in sectie 7.2.6 heeft deze parameter de beste trade-off tussen rekentijd en verbeterde waarde van de doelfunctie. Bij deze tests werd de waarde van deze parameter gevarieerd van 500 tot 5000 met tussenstappen van 500. Voor de rest werden de standaardwaarden uit tabel 7.2 gebruikt. De rekentijden voor het genetisch algoritme varieerden van 1 uur tot bijna 12 uur. Bij simulated annealing is α de logische keuze, zoals reeds werd beschreven in sectie 8.2.4. Opnieuw hebben alle andere parameterwaarden de staandaardwaarden zoals die in tabel 8.6
Hoofdstuk 9. Performantie van heuristieken
100
beschreven staan. De waarde van α werd gevarieerd van 99.5% tot 99.95% met tussenstappen van 0.05%. Dit resulteerde in berekeningstijden van 20 minuten tot meer dan 3 uur. 4.500.000
Simple heuristics performance
4.000.000
3.500.000
Objective function value
3.000.000
2.500.000
2.000.000
1.500.000
1.000.000
500.000
0 0
5.000
10.000
15.000
20.000
25.000
30.000
35.000
40.000
45.000
CPU Time (sec) Genetic Algorithm
Simulated Annealing
Figuur 9.2: Performantie van metaheuristieken voor verschillende rekentijden.
De resultaten van deze berekeningen worden getoond in grafiek 9.2. Het is snel duidelijk dat simulated annealing hier opnieuw beter gaat presteren dan het genetisch algoritme. Beide algoritmes bieden echter een zeer significante verbetering ten opzichte van de enkelvoudige heuristieken die als proxy voor de huidige planning worden beschouwd. De verbetering van de simulated annealing metaheuristiek vertoont in het begin een exponenti¨ele trend, wat aangeeft dat een licht hogere rekentijd aan te raden is wanneer deze voor handen is. De verbeteringen van het genetisch algoritme verloopt over heel het spectrum lineair. De voornaamste conclusie is hier dat het genetisch algoritme zelfs voor zeer sterk verhoogde rekentijd de oplossing van de simulated annealing heuristiek met standaardparameters niet kan verbeteren of evenaren. Grafiek 9.3 toont hoe de simulated annealing metaheuristiek naar een oplossing evolueert. Het is duidelijk dat voor een hogere waarde van α de zoekruimte veel breder is. Dit is zichtbaar door de grotere pieken bij de eerste iteraties, wat aangeeft dat de oplossing eerst veel verder
Hoofdstuk 9. Performantie van heuristieken
101
Solution code Sum of Value
10.000.000
9.000.000
8.000.000
Objective function value
7.000.000
6.000.000
5.000.000
Simple heuristics performance
4.000.000
3.000.000
2.000.000
1.000.000
0 0
1000
2000
3000
4000
5000
6000
7000
8000
9000
10000
11000
12000
13000
14000
15000
Iterations Solution number 99,50%
99,75%
99,95%
Dimension1
Figuur 9.3: Evolutie naar oplossing van simulated annealing heuristiek voor verschillende waarden van α
van de optimale oplossing evolueert alvorens ze opnieuw bij de optimale oplossing in de buurt gaat komen. Grafiek 9.3 toont dezelfde informatie voor het genetisch algoritme. De manier waarop dit algoritme naar een oplossing evolueert is veel gelijkmatiger dan de manier waarop het simulated annealing algoritme naar zijn oplossing convergeert. De verschillen tussen opeenvolgende iteraties zijn hier steeds groter dan bij het simulated annealing algortime, maar er is een veel grotere rekentijd nodig per iteratie. De gemiddelde verandering van de doelfunctie bedraagt voor het simulated annealing algoritme ongeveer 7.000 (α = 0.9995), voor het genetisch algoritme loopt dit op tot 25.000 (nakomelingen = 5000). Echter een enkele iteratie van het genetisch algoritme duurt 835 seconden, waar een iteratie van de simulated annealing heuristiek duurt gemiddeld minder dan 1 seconde duurt.
Hoofdstuk 9. Performantie van heuristieken
102
4.000.000
Simple heuristics performance
3.500.000
Objective function value
3.000.000
2.500.000
2.000.000
1.500.000
1.000.000
500.000
0 0
10
20
30
40
Iterations
500
2500
5000
Figuur 9.4: Evolutie naar oplossing van genetisch algoritme voor verschillende hoeveelheid nakomelingen per generatie.
9.3
Conclusie
Uit de bovenstaande analyse blijkt duidelijk dat het simulated annealing algoritme zowel efficienter als effectiever is om dit probleem op bedrijfsschaal te gaan optimaliseren. Op basis van de gegevens voorgesteld in grafiek 9.2 lijkt het optimaal om deze heuristiek een hogere alfa waarde te geven indien de rekentijd voorhanden is. De rekentijd zal bijgevolg toenemen maar is al bij al nog steeds vrij beperkt te houden. Op minder dan een uur kan vrijwel al het opwaarts potentieel worden bereikt voor problemen op bedrijfsschaal, dit is minder dan de tijd die nu dagelijks wordt gespendeerd aan het maken van de planning.
Hoofdstuk 10
Conflicten door beperkte beschikbaarheid van personeel Zoals eerder besproken in sectie 2.3 is het mogelijk dat omschakelingen op exact hetzelfde tijdstip worden gepland, wat kan zorgen voor een gebrek aan capaciteit voor de omschakeling. (Dit werd ge¨ıllustreerd in de onderste rij van figuur 2.2.) Vooral bij een grotere omvang van het probleem neemt de kans op dergelijke conflicten toe. Dit kan natuurlijk een nefast effect hebben op de doelfunctie als hierdoor jobs met vertraging worden opgeleverd. Deze dimensie werd nog niet opgenomen in de doelfunctie van de bovenstaande metaheuristieken, om de rekentijd te beperken. Het effect van deze conflicten wordt hier in detail bestudeerd. Dit hoofdstuk bestaat uit twee delen. Het eerste deel analyseert het voorkomen van conflicten bij geoptimaliseerde schema’s. Deze analyse gebeurt experimenteel op basis van 50 willekeurig gegenereerde datasets die van de huidige bedrijfssituatie starten. Het tweede deel van dit hoofdstuk behandelt de impact van het aantal machines op het voorkomen van dergelijke conflicten. Hier gebeurt de analyse niet langer op een experimentele manier maar wordt gebruik gemaakt van kansrekenen om voorspellingen te maken.
10.1
Analyse van conflicten bij geoptimaliseerde schema’s
Om te kunnen analyseren welke impact het beperkt aantal technici heeft werden 50 datasets gegenereerd die elk 300 jobs plannen op 32 machines die zich op twee productielocaties
103
Hoofdstuk 10. Conflicten door beperkte beschikbaarheid van personeel
104
bevinden. Voor deze datasets werd vervolgens een oplossing berekend met de genetische metaheuristiek uit hoofdstuk 7. In deze oplossingen werd vervolgens gezocht naar momenten waarop de vraag naar techniekers voor omschakelingen het aanbod overtreft. Hierbij werd verondersteld dat telkens een technieker stand-by is om omschakelingen uit te voeren. (Daar techniekers naast deze omschakelingen ook andere taken hebben is dit een realistische assumptie.) Concreet werd voor deze berekening gebruik gemaakt van tijdsincrementen van een kwartier. De eerste stap in dit proces was de volledige makespan (= de totale tijd nodig om alle orders in te plannen) van de oplossing verdelen in stukken ter grootte van dit increment. Voor al deze tijdsincrementen wordt bepaald hoe veel machines met een omschakeling bezig zijn. Dit wordt nagegaan door alle jobs op alle machines te overlopen, indien een omschakeling bezig is wordt een eenheid vraag toegevoegd aan het desbetreffende increment. Belangrijk hierbij is dat deze vraag wordt opgesplitst per productielocatie, daar omschakelend personeel natuurlijk eigen is aan de productielocatie. De output van deze berekening geeft het vraagpatroon weer per locatie. Door de deterministische productietijden laat dit toe te voorspellen wanneer de vraag het aanbod zal overtreffen. De resultaten tonen dat gegeven de huidige configuratie van de machines een gemiddeld schema conflicten ter waarde van 300 tijdsperiodes in de Belgische productielocaties en conflicten ter waarde van 860 tijdsperiodes in de Roemeense productielocaties zal veroorzaken. Het exacte moment waarop deze conflicten zich zullen voordoen is aan het begin van het schema met enige zekerheid te voorspellen, naar het eind van de productietijd toe zal deze voorspelling veel minder accuraat zijn. Dit heeft hier echter geen impact op de conclusies. Gemiddeld genomen zullen immers nog steeds het zelfde aantal conflicten voorvallen gedurende de productie. Ook voor de voorspelling van het moment waarop conflicten gaan plaatsvinden is dit niet noodzakelijk een groot probleem daar de planning niet voor de volledige periode wordt vastgelegd maar gedurende deze periode systematisch wordt herrekend en aangevuld met nieuwe orders. Dit betekent dat enkel de korte termijn voorspellingen van belang zijn voor het managen van de productiehal. Deze resultaten kunnen worden herrekend ( 1160 4∗24 ≈ 12) als een totale vertraging van ongeveer 12 dagen voor het volledige schema. (Ofwel: de gecumuleerde vertraging voor alle jobs bedraagt in totaal 12 dagen.) Dit betekent dat door het aanpassen van de personeelsplanning
Hoofdstuk 10. Conflicten door beperkte beschikbaarheid van personeel
105
om stilstand veroorzaakt door deze conflicten te vermijden netto 290 productie uren kunnen worden gewonnen in een planningshorizon van 6 weken.
10.2
Invloed van het aantal machines op het voorkomen van personeelsconflicten
Het voorkomen van personeelsconflicten is afhankelijk van het aantal machines dat zich op een bepaalde locatie bevindt. In deze sectie wordt berekend welk percentage van de tijd conflicten voorkomen gegeven het aantal machines. De eerste stap voor het berekenen van deze kans is de kans dat een individuele machine aan het omschakelen is. Deze kans (P (switch)) kan worden benaderd met behulp van de gemiddelde productietijd (Pavg ) en de gemiddelde tijd nodig voor een omschakeling (Savg ): P (switch) = Savg Pavg +Savg .
Wanneer de bovenstaande vergelijking wordt ingevuld gebruikmakend van echte
productiegegevens blijkt dat de waarde voor P (switch) hier gelijk is aan 2,25%. Dit betekent dat een machine gemiddeld genomen 2,25% van de bruto productietijd aan het omschakelen is. In de tweede stap van de berekening kan uitgegaan worden van een eenvoudige kansboom waarbij elke iteratie een machine voorstelt. Bij elke iteratie van de kansboom zijn er telkens twee opties: de volgende machine is aan het produceren of moet volgens de planning momenteel worden omgeschakeld. Deze opties hebben telkens de kansen zoals deze hierboven werden berekend. Hier wordt verondersteld dat onder het huidige management maximaal een omschakeling gelijktijdig plaatsvind omdat slechts een technieker beschikbaar is voor het uitvoeren van deze omschakelingen. Uitgaande van deze kansboom kan de kans worden berekend dat op een willekeurig moment meer dan 1 machine moet worden omgeschakeld volgens eender welke productieplanning die deze factor niet expliciet in rekening neemt. Deze kans wordt berekend in vergelijking 10.1, deze vergelijking maakt gebruik van de inverse om de berekening eenvoudiger te maken: P (conf lict) = 1 − P (no − conf lict) = 1 − [P (all − producing) + P (1 − switching)]. P (conf lict) = 1 − P (prod)m + m ∗ (P (switch) ∗ P (prod)m−1 )
(10.1)
In vergelijking 10.1 stelt m het aantal machines op de productielocatie voor, en P (prod) de
Hoofdstuk 10. Conflicten door beperkte beschikbaarheid van personeel
106
kans dat een machine aan het produceren is. Dit is de inverse kans van P (switch), bijgevolg is P (prod) = 1 − P (switch). Verder is P (prod)m de kans dat alle machines aan het produceren zijn en is m ∗ (P (switch) ∗ P (prod)m−1 de kans dat exact 1 machine aan het omschakelen is. Wanneer de assumptie dat slechts een machine gelijktijdig kan worden omgeschakeld wordt vervangen door de assumptie dat maximaal twee machines gelijktijdig kunnen worden omgeschakeld, kan de kans berekend wordenmet vergelijking10.2. Het verschil tussen vergelijking m−2 P (10.1) en vergelijking(10.2) is de term 1 + (m − n) ∗ (P (switch)2 ∗ P (prod)m−2 ) die de n=1
kans voorstelt dat exact twee machines op een willekeurig moment aan het omschakelen zijn. Indien de capaciteit voor omschakelen wordt verhoogd naar twee simultane omschakelingen zorgt deze situatie immers ook niet langer voor conflicten.
P (conf lict) = 1 − (P (prod)m + m ∗ (P (switch) ∗ P (prod)m−1 )+ " # m−2 X 1+ (m − n) ∗ (P (switch)2 ∗ P (prod)m−2 )) (10.2) n=1
70,00%
Chance of resource conflict
60,00%
50,00%
40,00%
30,00%
20,00%
10,00%
1 3 5 7 9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 45 47 49 51 53 55 57 59 61 63 65 67 69 71 73 75 77 79 81 83 85 87 89 91 93 95 97 99
0,00%
# Machines P(conflict)
P(dual conflict)
Figuur 10.1: Kans op tekort aan personeel door gelijktijdige omschakelingen op een willekeurig moment.
Wanneer deze kansen worden uitgerekend voor een verschillend aantal machines op een locatie
Hoofdstuk 10. Conflicten door beperkte beschikbaarheid van personeel
107
wordt grafiek 10.1 bekomen. Uit deze grafiek blijkt duidelijk dat naarmate de productiecapaciteit op een locatie wordt vergroot de kans op stilstand door tekort aan personeel sterk gaat toenemen. Voor de huidige situatie in Roemeni¨e kan worden gesteld dat gedurende 10% van de tijd machines onnodig staan te wachten op personeel voor een omschakeling. Met een logisch verlies aan productiecapaciteit tot gevolg. Het permanent verhogen van de capaciteit voor de omschakelingen lijkt voorlopig nog niet nodig. Het nodige personeel om twee simultane omschakelingen uit te voeren is immers beschikbaar, maar deze worden niet steeds prioritair toegekend aan omschakelingen. Het probleem kan momenteel worden opgelost door te anticiperen op eventuele omschakelingen en er dan ook voor te zorgen dat techniekers geen andere taken toegewezen krijgen wanneer zij nodig zijn voor het uitvoeren van omschakelingen. Echter wanneer deze productiefaciliteit verder in schaal zou toenemen, zal de impact van deze conflicten quasi lineair blijven toenemen tot op een moment dat het voordeliger wordt om naar de onderliggende curve over te stappen door het verhogen van de capaciteit.
10.3
Besluit
De bovenstaande analyse toont dat niettegenstaande het integreren van deze conflicten technisch niet haalbaar blijkt, het anticiperen van deze conflicten een merkbare verhoging van de productiecapaciteit kan opleveren. Concreet betekent dit dat een nieuwe planningsmodule naast het inplannen van de jobs ook steeds de conflicten betreffende omschakelingen zal gaan signaleren aan de productieleiding. De productieleiding is dan op zijn beurt verantwoordelijk om ervoor te zorgen dat de pieken in de vraag correct worden opgevangen en niet voor een verstoring van het productieproces gaan zorgen.
Hoofdstuk 11
Slotbeschouwing In dit werk werden verschillende methodes aangereikt om de productieplanning van een Belgische textielproducent te optimaliseren. Meer specifiek werden een exacte oplossingsmethode en twee heuristische oplossingsmethodes ontworpen en getest. Algemeen kan worden gesteld dat alle gebruikte methodes significant (tot 70%) meer economische waarde cre¨eren dan de huidige manuele methode - dewelke met enkelvoudige heuristieken werd benaderd. De exacte oplossingsmethode kan echter enkel worden toegepast op problemen van beperkte schaal. Voor problemen op bedrijfsschaal moet dus een keuze worden gemaakt tussen een van de twee ontwikkelde heuristieken. De simulated annealing metaheuristiek blijkt zowel effectiever (betere kwaliteit van de eindoplossing) als effici¨enter (sneller) te zijn dan het genetisch algoritme. Het verschil tussen deze technieken neemt eveneens verder toe naarmate de probleemgrootte toeneemt. Het integreren van personeelsconflicten gedurende de omschakelingen van de machines in de algoritmes bleek praktisch moeilijk haalbaar door de sterk verhoogde rekentijd die hiermee gepaard gaat. Het detecteren en anticiperen van deze conflicten kan echter wel resulteren in significante productiviteitswinsten (200 extra productie uren per maand, verdeeld over het volledige productieapparaat).
108
Deel IV
Bijlagen
109
Bijlage A
Uitbreiding op exacte oplossingsmethode: integratie van personeelsconflicten Deze bijlage beschrijft een uitbreiding op het model uit hoofdstuk 5 waarbij ook de vertragingen in de productie veroorzaakt door overschrijding van het aanbod van omschakelend personeel worden ge¨ıntegreerd. De structuur van dit hoofdstuk is analoog aan die van hoofdstuk 5. Opnieuw werd ook een extra hoofdstuk toegevoegd waarin de vergelijkingen verder worden ge¨ıllustreerd zodat hun werking volledig duidelijk is (Appendix B). Het hieronder beschreven model was echter te complex om te worden ge¨ımplementeerd. Meer duiding hierbij kan teruggevonden worden in sectie B.5.
110
Bijlage A. Uitbreiding op exacte oplossingsmethode: integratie van personeelsconflicten 111
A.1
Indices Index
Beschrijving
j
Job: j ∈ J = {0, ..., N J}
k
Positie op machine: k ∈ K = {0, ..., N K}
m
Machine: m ∈ M = {0, ..., N M }
t
Tijd (interval van 15 min): t ∈ T = {0, ..., N T }
p
Plaats (verschillende productielocaties): p ∈ P = {0, ..., N P }
Tabel A.1: Indices gebruikt in optimalisatiemodel
A.2 A.2.1
Variabelen Binaire variabelen Variable
Beschrijving
xjkm
= 1 als job j op positie k van machine m gepland is, anders 0.
lj
= 1 as job j te laat is, anders 0.
t fji
= 1 als omschakeling ji afgelopen is op moment t, anders 0.
btji
= 1 wanneer omschakeling ji kan starten, anders 0.
Tabel A.2: Binaire variabelen gebruikt in optimalisatiemodel
A.2.2
Integer variabelen
Variable
Beschrijving
tj
Aantal dagen dat job j te laat is.
t rji t ωpji ρtp
Personeel toegekend aan omschakeling ji op moment t. Personeel toegewezen aan omschakeling ji op moment t en plaats p. Aantal personeelsleden toegewezen aan omschakelingen op plaats p en moment t.
Tabel A.3: Integer variabelen gebruikt in optimalisatiemodel
Bijlage A. Uitbreiding op exacte oplossingsmethode: integratie van personeelsconflicten 112
A.2.3
Continue Variabelen Variable
Beschrijving
cj
Tijdstip waarop job j wordt opgeleverd.
λj
Vaste kosten gemaakt door laattijdig leveren van job j.
sji
Duurtijd van omschakeling ji.
Tabel A.4: Continue variabelen gebruikt in optimalisatiemodel
A.3
Parameters Parameter
Beschrijving
Fjm
Vaste kost wanneer job j te laat wordt opgeleverd op machine m.
Vj
Variabele kost per dag dat job j te laat wordt opgeleverd.
Pjm
Productietijd in dagen van job j op machine m.
Rj
Tijdstip waarop job j in het systeem wordt toegevoegd.
Dj
De deadline van job j in dagen.
Atp
Beschikbaar personeel voor omschakeling op tijdstip t en plaats p.
Wji
Aantal werkeenheden nodig voor een omschakeling van j naar i.
Lmp
Binaire waarde die de locatie van de machine aangeeft.
M
Een groot positief nummer.
NJ
Totaal aantal jobs dat ingepland moet worden.
NK
Maximaal aantal posities op een machine.
NM
Aantal machines.
NT
Bovengrens van beschouwde tijdsperiode.
NP
Aantal productielocaties.
Tabel A.5: Parameters gebruikt in optimalisatiemodel
noot: een werkeenheid wordt gedefinieerd als de hoeveelheid werk die een technieker gedurende een periode van 15 minuten kan presteren aan een omschakeling.
Bijlage A. Uitbreiding op exacte oplossingsmethode: integratie van personeelsconflicten 113
A.4
Doelfunctie M in.
NJ X
[(λj ) + (Vj ∗ tj )]
(A.1)
j=0
A.5
randvoorwaarden N M X NK X
j∈J
xjkm = 1,
(A.2)
m=0 k=0 NJ X
xjkm ≤ 1,
k ∈ K, m ∈ M
(A.3)
j=0 NJ X j=0
xjkm −
NJ X
xj(k−1)m ≤ 0,
k ∈ K\{1}, m ∈ M
(A.4)
j=0
cj ≥ Pjm ∗ xjkm ,
j ∈ J, m ∈ M, k = 1
ci + M (2 − xikm − xj(k−1)m ) ≥ Pim + sji + cj ,
cj ≥ Rj +
(A.5)
j ∈ J, i ∈ J, k ∈ K\{0}, m ∈ M (A.6)
N M X NK X
(Pjm ∗ xjkm ),
j∈J
(A.7)
m=0 k=0
(cj − Dj ) ≤ M lj ,
j∈J
tj ≥ (cj − Dj ) − (1 − lj )M,
λj ≥ Fjm ∗ lj − M (1 −
NK X
xjkm )
j∈J
j ∈ J, m ∈ M
(A.8)
(A.9)
(A.10)
k=0
Atp ≥ ρtp ,
t t ωpji ≥ Lmp ∗ rji − M ∗ (1 −
NK X
xjkm ),
t ∈ T, p ∈ P
t ∈ T, p ∈ P, j ∈ J, i ∈ J, m ∈ M
(A.11)
(A.12)
k=0
ρtp =
NJ X NJ X j=0 i=0
t ωpji ,
p ∈ P, t ∈ T
(A.13)
Bijlage A. Uitbreiding op exacte oplossingsmethode: integratie van personeelsconflicten 114
t rji ≤ 2,
NT X
t ∈ T, j ∈ J, i ∈ J
t rji ≥ Wji − M ∗ (2 − xikm − xj(k−1)m ),
(A.14)
j ∈ J, i ∈ J, k ∈ K\{0}, m ∈ M
(A.15)
t=0
NT X
t∗btji ≥ 24∗4∗cj −M ∗(2−xikm −xj(k−1)m ),
j ∈ J, i ∈ J, k ∈ K\{1}, m ∈ M (A.16)
t=0
NT X
t ∗ btji ≤ 24 ∗ 4 ∗ cj + 1 + M ∗ (2 − xikm − xj(k−1)m ),
j ∈ J, i ∈ J, k ∈ K\{1}, m ∈ M
t=0
(A.17) 0
t X t=0
0
t rji
≤
t X
j ∈ J, i ∈ J, t0 ∈ T
btji ∗ M,
(A.18)
t=0 NT X
btji ≤ 1,
j ∈ J, i ∈ J
(A.19)
t=0 0
t X
0
t rji
+ M (1 −
t X
t=0
t fji ) ≥ Wji ,
t0 ∈ T, j ∈ J, i ∈ J
(A.20)
t=0 0
0
t X t=0
t rji
t X t fji ) < Wji , − M(
t0 ∈ T, j ∈ J, i ∈ J
(A.21)
t=0 NT X
t fji ≤ 1,
j ∈ J, i ∈ J
(A.22)
t=0 NT NT X X t sji = (1/24) ∗ (1/4) ∗ [( (t ∗ fji )) + 1 − (t ∗ btji )], t=0
(A.23)
t=0
t xjkm , lj , fji , btji ∈ {1, 0},
t rji ∈ N,
j ∈ J, i ∈ J
j ∈ J, i ∈ J, k ∈ K, m ∈ M, t ∈ T
j ∈ J, i ∈ J, m ∈ M, t ∈ T
tj , cj , λj , sjim ≥ 0,
j ∈ J, i ∈ J, m ∈ M
(A.24)
(A.25)
(A.26)
Bijlage A. Uitbreiding op exacte oplossingsmethode: integratie van personeelsconflicten 115
A.6
Beknopte duiding bij vergelijkingen
Hier wordt een zeer korte omschrijving gegeven van de bovenstaande vergelijkingen. Voor een meer uitgebreide uitleg wordt verwezen naar hoofdstuk B. A.1 De doelfunctie minimaliseert de som van de variabele en de vaste kosten die worden gemaakt wanneer een job te laat wordt opgeleverd. De parameter λ bevat de vaste kosten, die tevens rekening houden met de locatie waarop geproduceerd wordt. De variabele kosten zijn vervat in het product van de variabele kost per dag en het aantal dagen dat te laat wordt opgeleverd (= de variabele tj ). A.2 Zorgt ervoor dat alle jobs die gepland moeten worden exact 1 positie op een van de machines worden toegewezen. A.3 Op een bepaalde positie van een machine mag maximaal 1 job worden gepland. (Voorkomt dat er meerdere jobs op eenzelfde positie van eenzelfde machine ingepland worden.) A.4 Er mogen geen gaten in de planning van een bepaalde machine zitten. Indien dit wel het geval is zorgt dit voor problemen voor het bepalen van de omschakelingen die plaatsvinden. Concreet betekent dit dat wanneer er een job gepland is op positie k van machine m, er tevens een job gepland moet zijn op positie k-1 van dezelfde machine. (Met uitzondering van de jobs die op de eerste positie van een machine ingepland staan.) A.5 De eindtijd van de job die op de eerste positie gepland staat is gelijk aan zijn productietijd op die machine. A.6 De eindtijd van de job is gelijk aan de eindtijd van de voorgaande job op dezelfde machine, plus de tijd nodig om de machine om te schakelen (inclusief de eventuele wachttijden op techniekers om de machine om te schakelen) en de tijd nodig voor het produceren van de job zelf. A.7 Een job kan niet worden geproduceerd alvorens deze in het systeem wordt opgenomen. (Release time constraint.) A.8 Waardetoekenning van de lateness parameter, deze moet gelijk zijn aan 1 wanneer een job te laat wordt opgeleverd en gelijk aan 0 wanneer dit niet het geval is. A.9 Waardetoekenning voor de tardiness parameter, die het aantal dagen dat een bepaalde
Bijlage A. Uitbreiding op exacte oplossingsmethode: integratie van personeelsconflicten 116 job te laat geleverd wordt weergeeft. (Deze informatie is nodig voor het berekenen van de variabele kosten.) A.10 Waardeberekening van de lambda variabele in de doelfunctie. De reden voor het introduceren van deze variabele is het lineariseren van deze doelfunctie: M in.
N PJ
[(
N M N K P P
[Fjm ∗ xjkm ] ∗ lj ) + (Vj ∗ tj )]
j=0 m=0 k=0
Het probleem met deze doelfunctie is dat het linkerdeel van de vergelijking een vermenigvuldiging van twee variabelen is en bijgevolg niet lineair is. Door het toevoegen van de extra variabele λ en de constraint A.10 wordt dit probleem opgelost. A.11, A.12 en A.13 Toegekende techniekers op een bepaald moment in de tijd mogen het beschikbare aantal techniekers op dat moment niet overschrijden. A.14 Er kunnen op een bepaald moment in de tijd niet meer dan 2 mensen helpen met het omschakelen van een machine. (1 persoon stelt de machine in, terwijl de andere persoon de bobijnen van de machine vervangt.) A.15 Indien de omschakeling doorgaat moeten exact voldoende werkuren worden toegewezen aan de omschakeling om deze te volbrengen. A.16 enA.17 De omschakeling jim kan slechts beginnen nadat de productie van job j volledig is afgerond. (= de completion time van job j.) A.18 Voordat een omschakeling begint (= voordat de productie van de vorige job is afgerond) mogen nog geen werkuren worden toegekend aan de omschakeling. A.19 De binaire variabele die die start aanduidt van een omschakeling kan voor een bepaalde omschakeling maximaal een maal gelijk worden aan 1. A.20 Deze vergelijking zorgt ervoor dat de waarde van de variabele die het einde van de t omschakeling aanduidt fjim nul is wanneer nog niet voldoende resources toegewezen zijn aan
de desbetreffende omschakeling. A.21 Hierin wordt afgedwongen dat de waarde van de variabele die het einde van de omschat keling aanduidt fjim 1 moet zijn op het moment dat er voldoende werkuren zijn toegewezen
aan de desbetreffende omschakeling.
Bijlage A. Uitbreiding op exacte oplossingsmethode: integratie van personeelsconflicten 117 A.22 De som van binaire variabele die het einde van de omschakeling aanduidt mag maximaal gelijk zijn aan 1 voor elke individuele omschakeling. A.23 Berekenen van de correcte waarde voor de omschakeltijd, inclusief de tijd die er moet gewacht worden op beschikbaar personeel. A.24 Geeft aan welke variabelen binair zijn. A.25 Geeft aan welke variabelen integer zijn. A.26 Geeft aan welke variabelen positief (maar niet binair) moeten zijn.
Bijlage B
Uitbreiding op exacte oplossingsmethode: Illustraties en duiding Dit hoofstuk geeft meer uitleg bij de werking van het model uit bijlage A. Het is vooral bedoeld om mensen die minder ervaring hebben met de interpretatie van MIP optimalisatiemodellen zicht te geven op de redeneringen achter de vergelijkingen. Dit hoofdstuk beschrijft eerst de verschillende indices, variabelen en parameters die worden gebruikt in het model. Daarna worden alle randvoorwaarden grondig besproken en waar nodig met een voorbeeld ge¨ıllustreerd. In dit hoofdstuk werden de indices, variabelen en vergelijkingen die reeds werden besproken in hoofdstuk 6 niet herhaald. Aan het eind van dit hoofdstuk worden de problemen bij de implementatie van de uitgebreide modelformulering in detail besproken.
B.1
Indices
De tijdsdimensie voor verschillende variabelen wordt weergegeven met behulp van de index t. Een increment van deze index staat gelijk aan een increment met 15 minuten. Daar de rest van de tijdstippen in het model in dagen worden gemeten zorgt dit voor een aantal die hieronder in meer detail worden beschreven. De maximale waarde van deze index moet kunnen garanderen dat zelfs de traagst mogelijke productieplanning beschouwd kan worden. 118
Bijlage B. Uitbreiding op exacte oplossingsmethode: Illustraties en duiding
119
Om dit te kunnen garanderen worden alle traagste productiesnelheden voor de verschillende jobs opgeteld, en vervolgens verhoogd met de traagste NJ omschakelingen, waarbij NJ het aantal jobs is dat in het model is opgenomen.
B.2
Beslissingsvariabelen
Hieronder worden de binaire, integer en continue variabelen uit het model kort besproken en indien nodig ge¨ıllustreerd.
B.2.1
De begin- en eindtijd variabelen
De binaire variabele btji geeft weer op welk moment de omschakeling van job j naar job i kan beginnen. Dit is het moment waarop de productie van job j volledig is afgelopen. Let wel dat dit niet het moment is waarop de omschakeling van de machine echt van start gaat, omdat de mogelijkheid bestaat dat er op dat moment geen personeel beschikbaar is om de omschakeling uit te voeren. Indien de variabele b44 3,6 gelijk is aan 1 betekent dit dat de omschakeling van job nummer 3 naar job nummer 6 kan plaatsvinden vanaf t = 44. Daar er gewerkt wordt met tijdsintervallen van 15 minuten betekent dit dat de omschakeling kan plaatsvinden 11 uur nadat de productie van start gaat. t geeft op analoge manier als de voorgaande variabele manier weer De binaire variabele fji
wanneer de omschakeling van job j naar job i op machine m afgelopen is. Dit is het moment waarop voldoende werkeenheden toegewezen zijn aan de omschakeling van de machine. (En deze dus volledig ingesteld is om het volgende product te beginnen produceren.)
Figuur B.1: Illustratie van het gebruik van de begin- en starttijd variabelen
Bijlage B. Uitbreiding op exacte oplossingsmethode: Illustraties en duiding
120
Figuur B.1 illustreert hoe uit de begin- en startvariabelen de werkelijke duurtijd van de omschakeling kan worden afgeleid. hoe dit concreet gebeurt wordt hieronder nogmaals ge¨ıllustreerd wanneer de randvoorwaarden worden besproken.
B.2.2
Werkblokken toekennen aan omschakeling
t geeft het aantal werkblokken weer dat wordt toegekend aan de omschaDe integer variabele rji
keling van job j naar job i op een bepaald moment t. Een werkblok wordt hierbij gedefinieerd als de hoeveelheid werk die een technicus in 15 minuten kan uitvoeren voor een omschakeling. Deze variabele wordt gebruikt om het moment te bepalen waarop er voldoende werk gebeurd is om de omschakeling als afgerond te beschouwen. Wanneer de som over t van de r variabelen voldoende is wordt de omschakeling als compleet beschouwd.
B.2.3
Exacte duurtijd van omschakeling
Om de exacte duurtijd van een omschakeling, inclusief de wachttijd op beschikbaar personeel, te bepalen wordt de variabele sji gebruikt. Deze variabele maakt samen met de eindtijd variabele cj de link tussen de t index die tijdsblokken van 15 minuten voorstelt en de andere tijdstippen van het model die (fracties van) dagen weergeven. De manier waarop de waarde van deze variabele wordt toegekend wordt hieronder uitgebreid besproken.
B.3
Parameters
Het gebruik van de verschillende parameters in het model wordt ge¨ıllustreerd aan de hand van een voorbeeld. Tabel ?? geeft de vaste kosten, variabele kosten, productietijden en deadlines voor tien jobs weer. Hoe deze waarden worden berekend kan teruggevonden worden in hoofdstuk 4. De tabel ?? toont dat afhankelijk van de machine die wordt gekozen om een bepaald product te produceren de vaste kosten en de productietijden verschillen. Wanneer de productietijd als waarde -1 heeft betekent dit dat de desbetreffende machine job j niet kan produceren. Dit betekent dat een andere machine geselecteerd zal moeten worden voor deze job. Doordat het opnemen van alle machines een probleem zou opleveren dat te moeilijk oplosbaar
Bijlage B. Uitbreiding op exacte oplossingsmethode: Illustraties en duiding
121
is worden in deze dataset slechts drie machines opgenomen. Dit heeft echter ook als gevolg dat niet alle producttypes kunnen worden opgenomen in het model, omdat het niet mogelijk is om een combinatie te maken van drie machines die alle producten kan produceren. In tabel B.1 wordt voor de tien jobs de matrix met de omschakeltijden weergegeven. Dit zijn de waarden die worden voorgesteld door parameter Wji . Parameter W39 krijgt in deze tabel de waarde 8 toegekend. Aangezien de index t met tijdsintervallen van 15 minuten werkt betekent dit dat een technicus twee uur (= 8 * 15 min) nodig heeft om deze omschakeling te vervolledigen. Hoe deze parameters gebruikt worden in het model zal worden ge¨ıllustreerd met behulp in sectie B.4.
B.4 B.4.1
Randvoorwaarden opgelegd aan het model Allocatie van schaarse grondstoffen
Hieronder worden de randvoorwaarden besproken die betrekking hebben tot het toewijzen van technici aan omschakelingen van machines. Dit is een problematiek die in de beschikbare literatuur over het plannen van parallel geschakelde machines veelal over het hoofd wordt gezien. Aan de hand van deze vergelijkingen wordt uiteindelijk de waarde van de sji variabele bepaald, die de duurtijd van de omschakeling weergeeft inclusief de tijd die machines moeten wachten op beschikbaar personeel voor de omschakeling. Hieronder worden randvoorwaarden A.11 tot en met A.23 besproken. randvoorwaarde A.11 tot en met A.13 Voor elke locatie waarop er geproduceerd wordt is er op elk moment in de tijd slechts een bepaald aantal mensen beschikbaar om de omschakelingen uit te voeren. De toegewezen hoeveelheden mogen aldus deze beschikbare hoeveelheden niet overschrijden voor een gegeven locatie. Dit wordt afgedwongen met behulp van randvoorwaarden A.11 tot en met A.13. Opnieuw wordt de werking van deze constructie ge¨ıllustreerd aan de hand van een voorbeeld. Veronderstel een situatie waarin jobs 1 en 9 gepland zijn op machine 1, jobs 2 en 5 gepland staan gepland staan op machine 2 en jobs 3 en 4 gepland staan op machine 3. Hierbij bevinden machine 1 en 2 zich op locatie 1 (L11 = 1 en L21 = 1) en machine 3 op locatie 2 (L32 = 1). Daarnaast wordt ook verondersteld dat al deze omschakelingen plaats kunnen vinden in het interval t = {1,...,10}.
Bijlage B. Uitbreiding op exacte oplossingsmethode: Illustraties en duiding
1
2
3
4
5
6
7
8
9
10
1
2
2
7
5
6
2
32
2
7
2
2
2
2
10
10
10
2
32
2
12
2
3
7
10
2
7
7
2
32
2
8
2
4
5
10
7
2
6
2
32
2
7
2
5
6
10
7
6
2
2
32
2
6
2
6
2
2
2
2
2
2
32
2
12
2
7
32
32
32
32
32
32
2
32
32
32
8
2
2
2
2
2
2
32
2
12
2
9
7
12
8
7
6
12
32
12
2
2
10
2
2
2
2
2
2
32
2
2
2
Tabel B.1: Omschakeltijden in aantal tijdsintervallen van 15 minuten.
t
At1
At2
t r19
t r25
t r34
1
2
2
2
0
2
2
2
2
2
0
2
3
2
2
2
0
2
4
2
2
1
2
1
5
2
2
0
2
0
6
2
2
0
2
0
7
2
2
0
2
0
8
2
2
0
2
0
9
2
2
0
0
0
10
2
2
0
0
0
Tabel B.2: Illustratie randvoorwaarde A.11
t.e.m. A.13
122
Bijlage B. Uitbreiding op exacte oplossingsmethode: Illustraties en duiding
123
De omschakelingen die plaats moeten vinden kunnen worden voorgesteld zoals in tabel B.2. In deze tabel is aan elke omschakeling voldoende personeel toegewezen om deze omschakeling volledig uit te voeren. Wanneer deze situatie echter in detail wordt bestudeerd is een conflict zichtbaar voor locatie 1 op t = 4, waar 3 eenheden personeel werden toegewezen terwijl er voor die locatie slechts 2 eenheden beschikbaar waren. Wanneer de vergelijkingen worden uitgeschreven wordt het onderstaande resultaat bekomen: Voor p = 1, t = 4
4 4 ω119 ≥ L11 ∗ r19 − M ∗ (1 −
NK X
4 4 x1k1 ) ⇒ ω119 ≥ 1 ∗ 1 − M ∗ (1 − 1) ⇒ ω119 ≥1
(B.1)
4 4 x2k2 ) ⇒ ω125 ≥ 1 ∗ 2 − M ∗ (1 − 1) ⇒ ω125 ≥2
(B.2)
k=1
4 ω125
≥ L21 ∗
4 r25
− M ∗ (1 −
NK X k=1
ρ41
=
NJ X NJ X
4 ω1ji ⇒ ρ41 = 3
(B.3)
j=1 i=1
A41 ≥ ρ41 ⇒ 2 ≥ 3
(B.4)
Vergelijkingen B.1 en B.2 stellen randvoorwaarde A.12 voor, ingevuld voor de omschakelingen die plaatsvinden op de eerste machine (dit zijn de omschakelingen van job 1 naar job 9 en job 2 naar job 5 respectievelijk). Het resultaat van deze vergelijkingen is dat de ω een bepaalde minimale waarde meekrijgt indien de omschakeling op een specifiek tijdstip en locatie plaatsvind. Indien dit niet het geval is kan de waarde van de ω variabele gelijk aan nul zijn. (Dit wordt echter niet expliciet afgedwongen daar het model enkel dient te controleren of er niet te veel personeel wordt toegewezen. Vergelijking B.3 gaat vervolgens deze informatie vervat in de ω variabelen omzetten naar de ρ variabele. Deze integer variabele aggregeert de waarde van ω over alle jobs. (Noot: de reden waarom hier een dergelijke constructie wordt gebruikt is opnieuw het bewaren van de lineariteit van de vergelijkingen.) Vervolgens wordt deze waarde van ρ gebruikt in randvoorwaarde A.11 die wordt uitgeschreven in vergelijking B.4. Hier wordt de randvoorwaarde niet gerespecteerd, daar drie personeelsleden worden toegewezen op locatie 1, waar er slechts 2 personeelsleden beschikbaar zijn op
Bijlage B. Uitbreiding op exacte oplossingsmethode: Illustraties en duiding
124
dat moment en die plaats. Doordat deze randvoorwaarde niet gerespecteerd wordt zal deze oplossing ook niet haalbaar zijn. Om de bovenstaande omschakelingen haalbaar te maken zou er een personeelslid minder moeten worden ingepland. Het haalbare schema voor deze omschakelingen wordt ge¨ıllustreerd in tabel B.3. randvoorwaarde A.14 Aan een bepaalde omschakeling kunnen maximaal twee personen gelijktijdig werkzaam zijn. Het toevoegen van een derde persoon heeft weinig zin en wordt in de praktijk nooit gedaan. Deze voorwaarde wordt afgedwongen door randvoorwaarde A.14 die stelt dat op eender welk moment voor eender welke omschakeling het aantal toegekende techniekers niet groter mag zijn dan twee. randvoorwaarde A.15 Wanneer een bepaalde omschakeling van job j naar job i doorgaat (afhankelijk van het productieschema) moeten voldoende werkuren worden toegewezen aan deze omschakeling. randvoorwaarde A.15 zorgt ervoor dat dit wordt afgedwongen. Veronderstel een situatie waarin twee jobs op een machine moeten worden geproduceerd, job 1 wordt op de eerste positie ingepland (x111 = 1) en job 2 wordt op de tweede positie ingepland (x221 = 1). Opnieuw wordt verondersteld dat het hier om jobs 1 en 2 uit tabel B.1 gaat, wat betekent dat aan de omschakeling van job 1 naar job 2 exact 2 werkeenheden moeten worden toegewezen. Dit zorgt voor de volgende invulling van randvoorwaarde A.15.
j = 1, i = 2, k = 2, m = 1
NT X
t r12 ≥ W12 − M ∗ (2 − x221 − x1(2−1)1 )
t=1
⇒
NT X t=1
t r12
≥ 2 − M ∗ (2 − 1 − 1) ⇒
NT X
t r12 ≥2−0
t=1
(B.5)
Bijlage B. Uitbreiding op exacte oplossingsmethode: Illustraties en duiding
125
t
At1
At2
t r19
t r25
t r34
t ω119
t ω219
t ω125
t ω225
t ω134
t ω234
ρt1
ρt2
1
2
2
2
0
2
2
0
0
0
0
2
2
2
2
2
2
2
0
2
2
0
0
0
0
2
2
2
3
2
2
2
0
2
2
0
0
0
0
2
2
2
4
2
2
1
1
1
1
0
1
0
0
1
2
1
5
2
2
0
2
0
0
0
2
0
0
0
2
0
6
2
2
0
2
0
0
0
2
0
0
0
2
0
7
2
2
0
2
0
0
0
2
0
0
0
2
0
8
2
2
0
2
0
0
0
2
0
0
0
2
0
9
2
2
0
1
0
0
0
1
0
0
0
1
0
10
2
2
0
0
0
0
0
0
0
0
0
0
0
Tabel B.3: Een correcte oplossing voor de illustratie van A.11 t.e.m. A.13
j = 2, i = 1, k = 2, m = 1
NT X
t r21 ≥ W21 − M ∗ (2 − x121 − x2(2−1)1 )
t=1
⇒
NT X
t r21 ≥ 2 − M ∗ (2 − 0 − 0) ⇒
t=1
NT X
t r21 ≥ −M
t=1
(B.6)
Vergelijking B.5 dwingt voor de omschakeling die plaats zal vinden af dat er voldoende personeel wordt toegewezen. Vergelijking B.6 toont dat randvoorwaarde A.15 geen hoeveelheid toegewezen werkeenheden zal afdwingen voor omschakelingen die niet plaatsvinden. randvoorwaarde A.16 en A.17 De variabele btji geeft het moment aan vanaf wanneer resources kunnen worden toegewezen aan een bepaalde omschakeling. Dit is het moment waarop de job op de voorgaande positie op de dezelfde machine is afgerond. randvoorwaarde A.16 en A.17 gaan deze binaire variabele op het correcte tijdstip t de waarde 1 meegeven. Doordat cj een continue variabele is die in aantal dagen wordt uitgedrukt is hier een vrij
Bijlage B. Uitbreiding op exacte oplossingsmethode: Illustraties en duiding
126
complexe constructie nodig om af te leiden op welk moment de omschakeling kan beginnen uitgedrukt in tijdsperiodes van 15 minuten (de index t). In deze vergelijking wordt cj vermenigvuldigd met 24 (# uren per dag) en vervolgens met 4 (# kwartieren per uur). Het getal dat wordt bekomen is echter nog steeds een kommagetal, waar er voor het toekennen van de variabele btji eigenlijk een natuurlijk getal nodig zou zijn. Dit wordt opgelost door twee vergelijkingen te gebruiken, die de waarde van de integer insluiten. Opnieuw wordt dit ge¨ıllustreerd aan de hand van een voorbeeld. Veronderstel dat job 1 en 2 uit tabel ?? achtereenvolgens op machine 1 worden gepland. (Respectievelijk op posities k = 1 en k = 2) Dit betekent dat de productie van job 1 afgerond is na 2.29 dagen. ( = de productietijd van order 1 op machine 1). Onderstaande vergelijkingen kunnen worden opgesteld:
j = 1, i = 2, k = 2, m = 1 NT X
t ∗ bt12 ≥ 24 ∗ 4 ∗ c1 − M ∗ (2 − x221 − x1(2−1)1 )
t=0 NT X
t ∗ bt12 ≥ 24 ∗ 4 ∗ 2.29 − M ∗ (2 − 1 − 1)
t=0 NT X
t ∗ bt12 ≥ 219.84
t=0
(B.7)
j = 1, i = 2, k = 2, m = 1 NT X
t ∗ bt12 ≤ 24 ∗ 4 ∗ c1 + 1 + M ∗ (2 − x221 − x1(2−1)1 )
t=0 NT X
t ∗ bt12 ≤ 24 ∗ 4 ∗ 2.29 + 1 + M ∗ (2 − 1 − 1)
t=0 NT X
t ∗ bt12 ≤ 220.84
t=0
(B.8)
Bijlage B. Uitbreiding op exacte oplossingsmethode: Illustraties en duiding
127
Vergelijking B.7 correspondeert met randvoorwaarde A.16, en vergelijking B.8 correspondeert met randvoorwaarde A.17. Het is duidelijk dat als gevolg van deze tangconstructie de waarde N PT van de som t ∗ bt12 gelijk moet zijn aan de waarde 220. Deze combinatie van randvoort=0
waarden zal er echter nog niet eenduidig voor zorgen dat dit dan ook b220 12 is die gelijk is aan 1, hiervoor wordt nog randvoorwaarde A.19 ge¨ıntroduceerd. Het gewenste resultaat van deze randvoorwaarde wordt nogmaals ge¨ıllustreerd in tabel B.4. Tevens kan worden opgemerkt dat opnieuw gebruik wordt gemaakt van de constructie M ∗ (2 − x221 − x1(2−1)1 ) die ervoor zorgt dat wanneer een omschakeling niet doorgaat de begin en eindtijden van deze omschakeling tussen -M en +M kunnen vallen, daar deze in dit geval geen nut hebben binnen het model. randvoorwaarde A.18 Een van de redenen waarom de bt12 variabele wordt gebruikt is voorkomen dat er personeel wordt toegekend aan een omschakeling alvorens deze werkelijk van start kan gaan. (Omdat het vorige order nog steeds in productie is.) Om dit te illustreren wordt opnieuw verder gebouwd op het voorbeeld dat hierboven werd gebruikt waarbij een omschakeling plaatsvind van job 1 naar job 2. Hierboven werd reeds de waarde van de bt12 variabele bepaald. De situatie in tabel B.5 is niet de correcte situatie zoals deze door het model wordt afgedwongen. Er wordt immers reeds personeel toegekend aan een omschakeling alvorens deze van start is gegaan. ( = alvorens het voorgaande order klaar is met produceren.) Dit geeft t
bt12
217
0
218
0
219
0
220
1
221
0
222
0
Tabel B.4: Illustratie randvoorwaarde A.16
Bijlage B. Uitbreiding op exacte oplossingsmethode: Illustraties en duiding t
t r12
bt12
217
0
0
218
2
0
219
0
0
220
2
1
221
2
0
222
1
0
223
0
0
224
0
0
225
0
0
226
0
0
128
Tabel B.5: Illustratie randvoorwaarde A.18
de onderstaande invullingen van randvoorwaarde A.18: j = 1, i = 2, t0 = 218 ⇒
218 X
t r12 ≤
t=1
0
j = 1, i = 2, t = 220 ⇒
220 X t=1
t r12
≤
218 X
bt12 ∗ M ⇒
t=1 220 X t=1
bt12
∗M ⇒
218 X
t r12 ≤0⇒2≤0
(B.9)
t r12 ≤1∗M ⇒4≤M
(B.10)
t=1 220 X t=1
Vergelijking B.9 toont aan dat de randvoorwaarde in dit geval wordt geschonden, de waarde van de som van de toegekende werkblokken is groter dan nul voordat de omschakeling werkelijk 218 gelijk aan nul moet worden. kan beginnen. Dit betekent concreet dat de waarde van r12
De tweede vergelijking B.10 daarentegen zorgt niet voor problemen, op dit moment mogen immers reeds werkuren worden toegekend aan de omschakeling. randvoorwaarde A.19 Deze randvoorwaarde gaat er in combinatie met randvoorwaarden A.16 en A.17 voor gaan zorgen dat de binaire variabele die het mogelijk begin van de omschakeling aanduidt de juiste waarde meekrijgt. In woorden stelt deze randvoorwaarde dat over de complete tijdsperiode de binaire variabele voor een specifieke omschakeling maximaal 1 keer gelijk aan 1 kan worden gesteld.
Bijlage B. Uitbreiding op exacte oplossingsmethode: Illustraties en duiding
129
Wanneer opnieuw de situatie wordt beschouwd waarin omgeschakeld wordt van job 1 naar job 2 op machine 1 kan indien randvoorwaarde A.19 niet ge¨ıntroduceerd wordt in het model zich een situatie voordoen waarin de btji variabele niet de correcte waarde meekrijgt. Concreet zou 120 dit bijvoorbeeld kunnen doordat zowel b100 12 en b12 gelijk worden gesteld aan 1. Dit zou geen
schending betekenen van randvoorwaarden A.16 en A.17 maar een foutief schema produceren. randvoorwaarde A.20 en A.21 t Deze randvoorwaarden controleren de allocatie van de correcte waarde aan de variabele fji
die het einde van de omschakeling van job j naar job i aanduid. Veronderstel opnieuw de situatie waarin van job 1 naar job 2 moet worden omgeschakeld op machine 1. De allocatie van werkblokken zoals in tabel B.6 zorgt ervoor dat de omschakeling afloopt op t = 221. Wanneer dit wordt ingevuld in de randvoorwaarden die de waarde voor de eindtijd variabele de correcte waarde meegeven wordt het onderstaande resultaat bekomen. (Noot: deze randvoorwaarde wordt steeds ingevuld voor j = 1 en i = 2.) t0 = 219 ⇒
219 X
t r12 + M (1 −
t=1
0
t = 219 ⇒
219 X
t f12 ) ≥ W12 ⇒ 0 + M (1 −
t=1 219 X t=1
t r12
219 X
t f12 )≥2
(B.11)
t=1
219 219 X X t t − M( f12 ) < W12 ⇒ 0 − M ( f12 )<2 t=1
(B.12)
t=1
De bovenstaande combinatie van randvoorwaarden A.20 en A.21 ingevuld voor t’ = 219 t variabele tot op dat moment gelijk is aan nul. (Afgedwongen verplicht dat de som van f12
door vergelijking B.11.) randvoorwaarde A.21 legt in dit geval geen strikte voorwaarde op, de waarde van de som mag voor deze randvoorwaarde zowel gelijk aan 0 als gelijk aan 1 zijn. t
t r12
bt12
219
0
0
220
1
1
221
1
0
222
0
0
Tabel B.6: Illustratie randvoorwaarde A.20 en A.21
Bijlage B. Uitbreiding op exacte oplossingsmethode: Illustraties en duiding
t0 = 220 ⇒
220 X
t r12 + M (1 −
t=1
220 X
t f12 ) ≥ W12 ⇒ 1 + M (1 −
t=1
t0 = 220 ⇒
220 X
t r12 − M(
220 X
t f12 )≥2
130
(B.13)
t=1 220 X
t=1
t f12 ) < W12 ⇒ 1 − M (
t=1
220 X
t f12 )<2
(B.14)
t=1
Wanneer t’ verhoogd wordt naar 220 wordt
220 P t=1
t gelijk aan 1, daar er reeds 1 werkeenheid r12
werd toegewezen aan deze omschakeling. Dit verandert voorlopig niets aan de interpretatie 220 P t van de bovenstaande vergelijkingen. Vergelijking B.13 dwingt af dat f12 gelijk is aan nul, t=1
wat nog steeds wordt toegelaten door vergelijking B.14. t0 = 221 ⇒
221 X
t r12 + M (1 −
t=1
221 X
t f12 ) ≥ W12 ⇒ 2 + M (1 −
t=1
t0 = 221 ⇒
221 X
t r12 − M(
221 X
t f12 )≥2
(B.15)
t=1 221 X
t=1
t f12 ) < W12 ⇒ 2 − M (
t=1
221 X
t f12 )<2
(B.16)
t=1
Nadat t’ nogmaals wordt ge¨ıncrementeerd zal randvoorwaarde A.20 (uitgeschreven in vergelij221 P t king B.15) er niet langer voor zorgen dat f12 gelijkgesteld wordt aan nul, deze vergelijking t=1
gaat nu immers op zowel voor waarde 1 als voor waarde 0 van deze uitdrukking. randvoorwaarde A.21 (uitgeschreven in vergelijking B.16) daarentegen zal nu gaan afdwingen dat 221 P t 221 de waarde 1 f12 de waarde 1 heeft. Concreet betekent dit dat de binaire variabele f12
t=1
krijgt. 0
t = 222 ⇒
222 X
t r12
+ M (1 −
t=1 0
t = 222 ⇒
222 X
t f12 )
≥ W12 ⇒ 2 + M (1 −
t=1 222 X t=1
t r12
− M(
222 X
t f12 )≥2
(B.17)
t=1 222 X t=1
t f12 )
< W12 ⇒ 2 − M (
222 X
t f12 )<2
(B.18)
t=1
Wanneer t’ verder wordt opgehoogd heeft dit geen effect op de interpretatie van de vergelijt>220 kingen. De waarde van alle opvolgende f12 variabelen zal gelijk aan nul worden gesteld.
randvoorwaarde A.22 Deze randvoorwaarde stelt dat de waarde van de variabele die het einde van de omschakeling aangeeft maximaal 1 keer over de gehele tijdsperiode de waarde 1 mag meekrijgen. De reden hiervoor is analoog aan de reden voor het introduceren van randvoorwaarde A.19.
Bijlage B. Uitbreiding op exacte oplossingsmethode: Illustraties en duiding
131
randvoorwaarde A.23 Deze randvoorwaarde bepaalt de tijd die nodig is voor het omschakelen van job j naar job i, inclusief de tijd die gewacht moet worden op beschikbaar personeel voor het uitvoeren van de omschakeling. Dit is het tweede doel van het begin en eindtijd variabelen. Om het gebruik van deze vergelijking te gaan illustreren wordt gebruik gemaakt van het voorgaande voorbeeld. (Het omschakelen van job 1 naar job 2 op machine 1 volgens de gegevens in tabel ??.) Dit voorbeeld wordt nogmaals samengevat in tabel B.7. Wanneer de randvoorwaarde wordt uitgeschreven voor j = 1 en i = 2 wordt het onderstaande resultaat bekomen.
NT NT X X t s12 = (1/24)∗(1/4)∗[( (t∗f12 ))+1− (t∗bt12 )] ⇒ s12 = (1/24)∗(1/4)∗[221+1−220] = 0.021 t=1
t=1
(B.19) De omschakeling van job 1 naar job 2 neemt in dit geval dus 0.21 dagen of omgerekend een half uur in beslag. De waarde van de sji variabele wordt vervolgens gebruikt in vergelijking A.6 om te bepalen wanneer order 2 kan afgerond worden.
B.4.2
Eigenschappen van variabelen bepalen
randvoorwaarden A.24 tot en met A.26 bepalen de limieten en de aard van de verschillende beslissingsvariabelen. Deze vergelijkingen zijn allemaal zeer eenvoudig en worden hier niet verder besproken. t
t r12
bt12
t f12
219
0
0
0
220
1
1
0
221
1
0
1
222
0
0
0
Tabel B.7: Illustratie randvoorwaarde A.23
Bijlage B. Uitbreiding op exacte oplossingsmethode: Illustraties en duiding
B.5
132
Problemen bij implementatie van optimalisatiemodel
Het model dat in dit en het voorgaande hoofstuk staat beschreven is in staat om de planningsproblematiek met volledige zekerheid optimaal op te lossen. Deze volledige zekerheid kan enkel gegarandeerd worden door een extreem groot aantal variabelen en beperkingen. Zelfs wanneer de infrastructuur voor handen zou zijn om een dergelijk model op te lossen zou dit nog steeds van weinig nut zijn voor het oplossen van de bedrijfsproblematiek. 10.000.000.000
1.000.000.000
# Variabelen
100.000.000
10.000.000
1.000.000
100.000
10.000
1.000
100
10
1 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50
# Jobs 2 Machines
3 Machines
4 Machines
30 Machines
Figuur B.2: Aantal variabelen nodig voor het oplossen van het exact optimalisatiemodel, uitgezet op logaritmische schaal
Figuur B.2 illustreert dat deze formulering zelfs bij zeer geringe omvang van het probleem een zeer groot aantal variabelen nodig heeft. Daarnaast toont de grafiek ook dat het aantal variabelen zeer snel toeneemt naarmate meer jobs ingepland moeten worden. Een interessante bemerking is ook dat de complexiteit slechts zwak toeneemt met het aantal machines.
Bibliografie D. Anghinolfi & M. Paolucci (2007). Parallel machine total tardiness scheduling with a new hybrid metaheuristic approach. Computers and Operations Research, 34(11):3471–3490. U. Bilge, F. Kirac, M. Kurtulan & P. Pekgun (2004). A tabu search algorithm for parallel machine total tardiness problem. Computers and Operations Research, 31(3):397–414. J. C. M. Carvalho & M. Ceccarelli (2001). A closed-form formulation for the inverse dynamics of a cassino parallel manipulator. Multibody System Dynamics, 5(2):185–210. G. Celano, A. Costa & S. Fichera (2008). Scheduling of unrelated parallel manufacturing cells with limited human resources. International Journal of Production Research, 46(2):405– 427. J.-F. Chen & T.-H. Wu (2006). Total tardiness minimization on unrelated parallel machine scheduling with auxiliary equipment constraints. Omega, 34(1):81–89. A. Gharehgozli, R. Tavakkoli-Moghaddam & N. Zaerpour (2009). A fuzzy-mixed-integer goal programming model for a parallel-machine scheduling problem with sequence-dependent setup times and release dates. Robotics and Computer-Integrated Manufacturing, 25(45):853–859. F. Glover & G. A. Kochenberger (2003). Handbook of metaheuristics. International series in operations research and management science. Kluwer Academic Publishers, Boston. 2002034082 edited by Fred Glover, Gary A. Kochenberger. ill. ; 25 cm. Includes bibliographical references and index. E. M. Goldratt & J. Cox (2004). The goal : a process of ongoing improvement. North River Press, Great Barrington, MA, 3rd rev. edition. 2006355019 by Eliyahu M. Goldratt and Jeff Cox. 23 cm. ”With interviews by David Whitford, editor at large, Fortune Small Business. Now with case study interviews--Cover. 133
Bibliografie
134
R. Graham, E. Lawler, J. Lenstra & A. Kan (1979). Optimization and approximation in deterministic sequencing and scheduling: a survey. 5:287 – 326. ISSN 0167-5060. URL http://www.sciencedirect.com/science/article/pii/S016750600870356X. A. Grigoriev, M. Sviridenko & M. Uetz (2007). Machine scheduling with resource dependent processing times. Mathematical Programming, 110(1):209–228. D.-W. Kim, K.-H. Kim, W. Jang & F. Frank Chen (2002). Unrelated parallel machine scheduling with setup times using simulated annealing. Robotics and Computer-Integrated Manufacturing, 18(3-4):223–231. K. Kumar, V. Selladurai, K. Raja & K. Elangovan (2011). Ant colony approach for makespan minimization on unrelated parallel machines. International Journal of Engineering Science and Technology, 3(4):3113–3120. Kurz & Askin (2001). Heuristic scheduling of parallel machines with sequence-dependent set-up times. International Journal of Production Research, 39(16):3747–3769. B. Mirkin & I. Muchnik (2009). Handbook of combinatorial optimization. International Journal of Production Research. M. Pfund, J. W. Fowler & J. N. D. Gupta (2008). A survey of algorithms for single and multi-objective unrelated parallel-machine deterministic scheduling problems. N. Piersma & W. van Dijk (1996). A local search heuristic for unrelated parallel machine scheduling with efficient neighborhood search. Mathematical and Computer Modelling, 24(9):11–19. R. Ruiz & A. Allahverdi (2007). No-wait flowshop with separate setup times to minimize maximum lateness. International Journal of Advanced Manufacturing Technology, 35(5):551– 565. M. Sewell, J. Samarab, R. Rodrigo & K. Mcisaac (2008). The rank-scaled mutation rate for genetic algorithms. R. Tavakkoli-Moghaddam, F. Taheri, M. Bazzazi, M. Izadi & F. Sassani (2009). Design of a genetic algorithm for bi-objective unrelated parallel machines scheduling with sequencedependent setup times and precedence constraints. Computers amp; Operations Research, 36(12):3224–3230.
Bibliografie
135
E. Vallada & R. Ruiz (2011). A genetic algorithm for the unrelated parallel machine scheduling problem with sequence dependent setup times. European Journal of Operational Research, 211(3):612–622. M. Yazdani, M. Amiri & M. Zandieh (2010). Flexible job-shop scheduling with parallel variable neighborhood search algorithm. Expert Systems with Applications, 37(1):678–687.