# Simplexmethode

De simplexmethode (of het simplexalgoritme) is een methode in de wiskundige optimalisatie (zie ook operationeel onderzoek). De techniek werd in 1947 door George Dantzig ontwikkeld. De simplexmethode lost een lineaire-optimaliseringsprobleem in een eindig aantal stappen op, of stelt de onoplosbaarheid van het probleem vast. In theoretische gevallen kunnen cycli optreden, die het vinden van de optimale oplossing verhinderen. De naam komt van het feit dat de vergelijkingen van het probleem een simplex beschrijven, waarvan de rand bij het vinden van de oplossing beschreven wordt.

## Geschiedenis

De Russische wiskundige Leonid Kantorovich behandelde in 1939 als eerste lineaire optimalisatie in zijn boek "Wiskundige Methoden in de Organisatie en Planning van de Productie". Kort daarna presenteerde de Amerikaan F.L. Hitchcock een werk over een transportprobleem. Toen erkende men echter nog niet de betekenis van dit werk.
De doorbraak voor de simplexmethode kwam er door George Dantzig in 1947, toen ook hij een werk over lineaire optimalisatie publiceerde. In eerste instantie toonden het Amerikaanse Leger, in het bijzonder de US Air Force, interesse in dit werk - men wilde de militaire acties optimaliseren. De daaropvolgende jaren ontwikkelden John von Neumann en Oskar Morgenstern het algoritme verder.
Met de opkomst van computers halverwege de jaren 50 kon men grotere problemen oplossen. Men ontwikkelde speciale varianten van de simplexmethode die op de eigenschappen van de computers uit die tijd afgesteld waren. Ze gingen bijvoorbeeld heel zuinig om met het gebruik van het hoofdgeheugen.
Gelijktijdig - na 1950 - ontdekte de industrie dat men met de simplexmethode ook optimalisatieproblemen in hun domein kon oplossen. In het bijzonder gebruiken olieraffinaderijen de methode. In 1960 optimaliseerde W. Knödel de toelevering van de suikerfabrieken in Oostenrijk. Daarmee verlaagde hij de transportkosten met ongeveer 10 procent.

## Definitie van het lineair optimalisatieprobleem

    Maximaliseer de doelfunctie $G:\mathbb {R} ^{n}\to \mathbb {R} ,\;x\mapsto a^{T}x$
    onder de nevenvoorwaarden $Ax\leq b$ en $x\geq 0$.
       Daarbij zijn
          A: m×n-matrix
          a: n-dimensionale vector
          b: m-dimensionale vector
       gegeven parameter, en de doelfunctie
       $a^{T}x=\langle a,x\rangle =\sum _{i=1}^{n}a_{i}x_{i}$ is als scalair product een lineaire functie in $x$.
De doelfunctie wordt ook wel objectief genoemd, de nevenvoorwaarden ook wel restricties.
$Ax$ is hier de notatie van een matrixvermenigvuldiging. De 'ongelijkheden' $Ax\leq b$ en $x\geq 0$ moeten component per component opgevat worden, bijvoorbeeld $x\geq 0$, met $x=(x_{1},x_{2},\dots ,x_{n})^{T}$, betekent hetzelfde als $x_{1}\geq 0,\dots ,x_{n}\geq 0$.
Die punten (vectoren) $x\,$ die aan alle nevenvoorwaarden tegelijkertijd voldoen (dit zijn geldige oplossingen van het probleem) vormen een convexe polyeder in de n-dimensionale ruimte, een zogeheten simplex. Wanneer er geen geldige oplossingen zijn, spreken de nevenvoorwaarden elkaar tegen en is het probleem niet oplosbaar.
Wanneer de verzameling geldige oplossingen niet leeg is, zijn er twee mogelijkheden:
* De doelfunctie G heeft in de toegelaten punten $x\,$ een bovengrens: Hier kan men zien dat het maximum van de doelfunctie door een punt op de rand van de geldige verzameling ingenomen wordt, en wel door een hoekpunt. Het simplexalgoritme verplaatst zich, vertrekkend van een geldige beginoplossing, langs de hoekpunten van de polyeder naar zijn optimale oplossing.
* De doelfunctie G is niet naar boven begrensd: In dat geval stelt men vast dat er een oneindig lange zijde van de polyeder met oneindige toenemende doelfunctie bestaat. Men zegt dat de optimale waarde van de doelfunctie $+\infty$ is, en men kan geen optimale oplossing daarvoor geven.
Wanneer geen geldige beginoplossing gekend is, kan men een lichtjes aangepast lineair optimaliseringsprobleem gebruiken, om een beginoplossing te vinden, en voor de vaststelling van de oplosbaarheid van het probleem. (De zogenoemde Fase I met kunstmatige variabelen).
Bij veel probleemstellingen komt het voor dat men de doelfunctie G niet wil maximaliseren, maar wil minimaliseren. Bijvoorbeeld moeten bij een transportprobleem de transportkosten zo klein mogelijk zijn, of moet de wegen zo kort mogelijk zijn. In dit geval voert men eenvoudigweg een nieuwe doelfunctie $H=-G$ in, en men maximaliseert deze.

## Voorbeeld

Aan de hand van een eenvoudig voorbeeld wordt de oplossingsmethode stap voor stap getoond. Reële problemen kunnen meer dan duizend randvoorwaarden hebben, waar men niet onmiddellijk eventuele tegenstellingen in kan herkennen, of dus niet direct het bestaan van een geldige beginoplossing kan vinden.

### Probleemstelling in woorden

Een firma biedt 2 verschillende producten aan. Er staan 3 machines A, B en C ter beschikking. Machine A heeft een maximale maandelijkse looptijd (capaciteit) van 170 uur, machine B van 150 uur en machine C van 180 uur. Een eenheid van product 1 geeft een winstbijdrage van 300 euro, een eenheid van product 2 daarentegen 500 euro. De vaste kosten bedragen 36.000 euro per maand. Als men 1 eenheid van product 1 vervaardigt, dan gebruikt men daarvoor 1 uur machine A en daarna 1 uur machine B. 1 eenheid van product 2 vereist achtereenvolgens 2 uur machine A, 1 uur machine B en 3 uur machine C.

### Wiskundige formulering met ongelijkheden

Stel dat het bedrijf per maand x$_1$ eenheden van product 1 en x$_2$ eenheden van product 2 vervaardigt, dan bedraagt de opbrengst (winst of verlies)
Deze opbrengst wil het bedrijf maximaliseren.
Aangezien de capaciteit van de machines begrensd is, kan men geen willekeurig grote hoeveelheden produceren. Uit de bezettingstijden leidt men de volgende nevenvoorwaarden voor x$_1$ en x$_2$ af:
Negatieve productiehoeveelheden zijn niet mogelijk, dus gelden nog de nevenvoorwaarden

### Omzetten in gelijkheden

Omdat het met stelsels gelijkheden gemakkelijker rekenen is, worden de ongelijkheden in gelijkheden omgezet. Daartoe voert men de zogenoemde slackvariabelen y$_A$, y$_B$ en y$_C$ in, die de niet benutte tijd van de respectievelijke machines voorstellen. Men vormt bijvoorbeeld de nevenvoorwaarde voor machine A in volgende gelijkheid om
Vanzelfsprekend mogen de slackvariabelen niet negatief zijn (niet-negativiteitsvoorwaarde).
Het probleem kan nu op volgende manier geformuleerd worden:
Maximaliseer de doelfunctie G onder de volgende nevenvoorwaarden:

### Het Simplex-Tableau

Deze vergelijkingen plaatst men in een zogeheten simplex-tableau, dat een verkorte voorstelling van het bijhorende stelsel van vergelijkingen is. Alle volgende bewerkingen zijn wiskundige berekeningen met het stelsel vergelijkingen.
De variabelen in de titelrij (x$_1$, x$_2$) heten niet-basisvariabelen, de variabelen in de eerste kolom basisvariabelen (Y$_A$, Y$_B$, Y$_C$). De getallen in "rij G" - de vergelijking van de doelfunctie - zijn de coëfficiënten van de doelfunctie. De variabelen b$_1$, b$_2$ en b$_3$ bevatten de waarde van de rechterleden.

### Bepaling van een geldige beginoplossing

Vervolgens bepaalt men een geldige oplossing. Dat zijn waarden voor x$_1$ en x$_2$ die aan alle nevenvoorwaarden voldoen.
Men kan in dit geval onmiddellijk een triviale oplossing geven, namelijk
Deze oplossing wordt in het bovenstaande tableau opgenomen. x$_1$ en x$_2$ zijn de niet-basisvariabelen en hebben steeds de waarde "0". De basisvariabelen hebben de waarde van de onbenutte machinelooptijden. Aangezien met de beginoplossing niets geproduceerd wordt: Y$_A$ = 170, Y$_B$ = 150, Y$_C$ = 180. Deze beginoplossing is natuurlijk onbevredigend: er wordt niet geproduceerd, de opbrengst G bedraagt -36.000 euro, wat verlies betekent. Daarom moet men proberen om betere geldige oplossingen te vinden.

### Verbetering van de oplossing door middel van een simplex-iteratie

In een iteratie van het simplexalgoritme wordt een basisvariabele met een niet-basisvariabele omgewisseld. Dit is een verwisselingsstap.
In aanmerking komen niet-basisvariabelen die een negatieve coëfficiënt in de doelfunctie hebben. Tussen deze variabelen zoekt men die variabele die bij omwisseling de grootste toename van de doelfunctie levert. Stel k het nummer van de kolom van de te verwisselen niet-basisvariabele, en r het rijnummer van de te verwisselen basisvariabele. De verwisselingsstap komt overeen met juist één stap bij het oplossen van een stelsel van gelijkheden, waarbij men rij r oplost naar de variabelen x$_k$ en dan x$_s$ in de rest van de gelijkheden inzet. Stel a$_{ij}$ een matrix-element uit de simplex-tableau, dan heet a$_{rk}$ het Pivotelement van de tableau. Kolom k is de pivotkolom, rij r de pivotrij.
In de verwisselingsstap berekent men het nieuwe tableau op volgende wijze:
Pivotelement:
       $a_{rk}:={1 \over a_{rk}}$ (Formule 1)
Pivotrij voor j $\neq$ k:
       $a_{rj}:={a_{rj} \over a_{rk}}$ (2)
       $b_{r}:={b_{r} \over a_{rk}}$ (3)
Pivotkolom voor i $\neq$ r:
       $a_{ik}:={-{a_{ik} \over a_{rk}}}$ (4)
overige elementen:
       $a_{ij}:={a_{ij}-{a_{ik}a_{rj} \over a_{rk}}}$ (5)
       $b_{j}:={b_{j}-{a_{ik}a_{rj} \over a_{rk}}}$ (6)

#### Keuze van de pivotelementen

Formule (6) is bepalend voor de keuze van de pivotrij r en de pivotkolom k. Men zoekt die waarden voor r en k waarvoor de coëfficiënt van de doelfunctie met index k negatief is, en waarvan de waarde van
    $G:=G-{a_{0k}\,b_{j} \over a_{rk}}$
zo groot mogelijk wordt. Daarbij dient men op te merken dat na de verwisseling de ontstane oplossing nog steeds geldig moet zijn, dit betekent dat alle nevenvoorwaarden vervuld moeten blijven. Daartoe kiest men r en k zodanig dat het quotiënt
    ${b_{j} \over a_{rk}}$ met $a_{rs}\neq 0$
binnen een kolom minimaal wordt.
Kolom 1 heeft een negatieve coëfficiënt van de doelfunctie van -300, en komt als pivotkolom k in aanmerking. Berekening van de coëfficiënten van deze kolom 1:
Het kleinste quotiënt 150 krijgt men dus in rij 2. Met pivotelement a$_{2,1}$ berekent men de nieuwe waarde van de doelfunctie als
    $G=-36000-{(-300)\cdot 150 \over 1}=9000$.
Ook kolom 2 heeft met -500 een negatieve coëfficiënt in de doelfunctie, en komt eveneens als pivotkolom k in aanmerking. De berekeningen van de coëfficiënten van deze kolom 2:
Het kleinste quotiënt 60 wordt gevonden in rij 3. Met het mogelijke pivotelement a$_{3,2}$ berekent men de nieuwe waarde van de doelfunctie als
    $G=-36000-{(-500)\cdot 180 \over 3}=-6000$.
Het is dus gunstiger om als pivotelement het eerste getal uit rij r = 2 en kolom k = 1 te kiezen, aangezien men daarmee de doelfunctie het meest doet toenemen, namelijk van -36000 naar 9000.

#### Uitvoeren van een verwisselingsstap

In de verwisselingsstap wordt de basisvariabele Y$_B$ met de niet-basisvariabele x$_1$ verwisseld.
Het simplex-tableau wordt volgens de bovenstaande regels omgerekend.

#### Het nieuwe simplex-tableau

Na omrekening verkrijgt men een nieuw simplex-tableau:
De variabele x$_1$ is tot de basis toegetreden, de variabele Y$_B$ verschijnt in de titelrij en is dus een niet-basisvariabele.
Deze oplossing betekent: er worden 150 eenheden van product 1 gemaakt (rechterlid van de rij met x$_1$). Van product 2 wordt niets vervaardigd (x$_2$ is een niet-basisvariabele). Daarmee heeft het bedrijf een opbrengst van 9000 euro. Machine A staat 20 uur per maand stil (draait slechts 150 van de 170 mogelijke uren). Machine B draait op volle toeren. Machine C staat 180 uur stil, ze wordt dus helemaal niet gebruikt. Dit kan men ook afleiden uit de opgave: machine C wordt enkel bij de productie van product 2 gebruikt. Aangezien dit product niet vervaardigd wordt, gebruikt men machine C dan ook niet.

### Verdere verbetering van de oplossing

Aangezien de doelfunctie in het nieuwe simplex-tableau nog een negatieve coëfficiënt heeft, kan men die oplossing verder verbeteren. Dit gebeurt door een volgende simplex-iteratie.
Bij de keuze van de pivotelementen komt enkel kolom k = 2 in aanmerking, aangezien enkel hier de coëfficiënt van de doelfunctie negatief is. De pivotrij kiest men volgens het minimum van de coëfficiënten "rechterlid/kolom 2", dus het minimum uit
Men kiest dus rij r = 1 als pivotrij. Het pivotelement is a$_{1,2}$ = 1, en de basisvariabele Y$_A$ wordt met de niet-basisvariabele x$_2$ verwisseld.
Na analoge berekeningen als hierboven wordt het nieuwe simplex-tableau:

### De optimale oplossing

Aangezien de doelfunctie geen negatieve coëfficiënten meer bevat, is de optimale oplossing bereikt.
Er worden 130 eenheden van product 1 en 20 eenheden van product 2 gefabriceerd. Daarmee behaalt het bedrijf een winst van 13.000 euro. Machine A en B worden ten volle benut. Machine C draait 60 uur, en heeft dus nog een onbenutte capaciteit van 120 uur (in het tableau is Y$_C$ = 120).

### Opmerking over de keuze van het pivotelement

In het voorbeeld wordt het pivotelement zodanig gekozen, dat de waarde van de doelfunctie bij de iteratiestap zo groot mogelijk wordt. Dit noemt men de Greatest-Change-Methode.
Alternatief kan men ook die kolom als pivotkolom nemen, waarvoor de coëfficiënt van de doelfunctie negatief en in absolute waarde het grootst is. Aangezien deze coëfficiënt een soort toename voorstelt, spreekt men van de Steepest-Unit-Ascent-Methode.
Vaak - maar niet altijd - levert de Greatest-Change-Methode met weinig simplex-iteraties een optimale oplossing.

### Grafische oplossing

In dit eenvoudige voorbeeld kent het probleem twee variabele parameters, in dat geval kan de oplossing grafisch voorgesteld worden in een 2-dimensionale figuur. Op de figuur worden de ongelijkheden voorgesteld als halfvlakken. De grensrechte wordt telkens afgebeeld. De doorsnede van deze halfvlakken bepaalt de regio met geldige oplossingen. Dit gebied is op de figuur geel gekleurd. De doelfunctie vormt ook steeds een rechte waarop het objectief eenzelfde waarde heeft, op de figuur zijn de rechten weergeven voor G = 0 en G = 13.000, in feite kan men een continue verzameling rechten tekenen. Bemerk dat al deze rechten die de doelfunctie voorstellen evenwijdig zijn, hoe groter G, hoe verder de rechte van de oorsprong ligt in dit geval. Men kan dan ook besluiten dat men zal wil nastreven deze rechte zo ver mogelijk van de oorsprong te krijgen, maar de rechte moet wel nog een punt binnen de geldige regio hebben, aangezien aan de nevenvoorwaarden moet voldaan zijn. Aangezien deze regio een convexe figuur is, ligt het maximum op een van de hoekpunten, dit is het snijpunt met de rechte waar G = 13.000, waar men x$_1$ = 130 en x$_2$ is 20 kan aflezen.
Op de figuur valt ook de werking van de simplex-iteraties grafisch te volgen. Het algoritme begon met beginoplossing (x$_1$=0, x$_2$=0), wat overeenkomt met het punt A. Dit is tevens het snijpunt van de rechten bepaald door de niet-negativiteitsvoorwaarden, namelijk de rechten met vergelijkingen x$_1$=0 en x$_2$=0. In tweede iteratie werd x$_1$ aan de basis toegevoegd, Y$_B$ werd nu een niet-basisvariabele. Wanneer men nu het snijpunt neemt van de rechten behorend bij deze niet-basisvariabelen, namelijk de rechte x$_2$ = 0 en de rechte B, dan bekomt men het punt B(150,0), en dit blijkt inderdaad ook de oplossing te zijn die men vindt bij de tweede iteratiestap. Bemerk dat vanuit punt B, men de doelfunctie enkel kan doen toenemen wanneer men zich verplaatst langs rechte B, wanneer men zich verplaatst langs de x-as, daalt de doelfunctie. Het is dan ook inderdaad Y$_B$ die uit de basis gehaald werd. In de derde iteratiestap zijn dan ook Y$_A$ en Y$_B$ niet-basisvariabelen: het snijpunt van de bijhorende rechten is opnieuw de oplossing van deze iteratiestap. Op deze manier loopt het simplex-algoritme langs de zijden van het geldig gebied, zolang er een zijde beschikbaar is langswaar de doelfunctie stijgend is, tot een maximum bereikt is.

## Alternatief tableau

Het simplex-tableau kan op een alternatieve manier ingedeeld worden, waarbij men alle variabelen (zowel basis- als niet-basisvariabelen) in de hoofding plaatst. Daaronder noteert men dan de doelfunctie en de nevenvoorwaarden, komt een variabele niet voor, dan schrijft men er een 0.
Het oorspronkelijke tableau ziet er dan als volgt uit:
We zien dat de rijen/kolommen Y$_{a,b,c}$ een eenheidsmatrix vormen. Deze tabel toont de oplossing { x$_1$ = 0, x$_2$ = 0, Y$_A$ = 170, Y$_B$ = 150, Y$_C$ = 180 } . Analoog aan voorheen kiest men kolom 2, omdat -500 in absolute waarde het grootst is, en een wijziging van x2 dus een potentieel grotere verandering in de objectieffunctie kan veroorzaken. Men kiest als pivotrij, rij YC. Dit besluit men opnieuw uit de berekening van rechterlid/x$_2$. Voor rij YA bekomt men 170/2 = 85, voor rij YB 150/1 = 150, voor rij YC 180/3 = 60. Men ziet echter dat wanneer men voor x$_2$ de waarde 85 zou gebruiken, de nevenvoorwaarde van rij YC nooit meer voldaan zou kunnen worden. Immers, vult men x = 85 in in 3 x$_2$ + Y$_C$ = 255 + Y$_C$, dan ziet men dat deze som nooit kleiner kan worden dan 180, aangezien Y$_C$ positief is, en de ongelijkheid zou dus nooit voldaan kunnen worden. Daarom kiest men die rij waar de verhouding rechterlid/x$_2$ het kleinst is, of dus rij YC. Ook hier voert men nu elementaire rijbewerkingen op de matrix uit, wat leidt tot het volgende tableau:
Opnieuw is te zien dat x$_2$ nu een basisvariabele is. De eenheidsmatrix is ook weer terug te vinden, met Y$_{a,b}$ en x$_2$ In de voorlopige oplossing heeft x$_2$ de waarde 60. De eerste kolom heeft nu een negatief element, rij YA heeft de kleinste kolomverhouding, pivoteren levert:
Aangezien er nog één kolom met negatieve coëfficiënt blijft, levert een laatste maal pivoteren met rij YB:
Men leest uiteindelijk dezelfde oplossing x$_1$ = 130 en x$_2$ = 20 af.

## Algoritme in pseudocode

In deze code wordt de begintabel als de matrix "tab" voorgesteld, waarbij G de onderste rij vormt en niet de bovenste. De matrix ziet er voor het voorbeeld als volgt uit: