1 Szakdolgozat Miskolci Egyetem Nemlineáris programozás Készítette: Horváth Gábor Programtervező informatikus hallgató Témavezető: Dr. Nagy Tamás egye...
Készítette: Horváth Gábor Programtervező informatikus hallgató Témavezető: Dr. Nagy Tamás egyetemi docens, Alkalmazott Matematikai Tanszék
Miskolc, 2013
Miskolci Egyetem Gépészmérnöki és Informatikai Kar Alkalmazott Matematikai Tanszék
Szám:
Szakdolgozat Feladat Horváth Gábor (GWUVPQ) programtervező informatikus jelölt részére. A szakdolgozat tárgyköre: operációkutatás A szakdolgozat címe: Nemlineáris programozás A feladat részletezése: Irodalomkutatás. Történeti áttekintés. Feltételes és feltétel nélküli optimalizálás. Newton és kvázi Newton módszerek, SUMT módszer. Bemenő adatok beviteli terve, kimenő adatok megjelenítési terve. Az algoritmusok programkódjának elkészítése. WEB-es megjelenítés.
Témavezető: Dr. Nagy Tamás, egyetemi docens A feladat kiadásának ideje: 2011.
................................. szakfelelős 2
1. szükséges (módosítás külön lapon) A szakdolgozat feladat módosítása nem szükséges ......................
...........................
dátum
témavezető(k)
2. A feladat kidolgozását ellenőriztem: témavezető (dátum, aláírás):
konzulens (dátum, aláírás):
..............
.............
..............
.............
..............
.............
3. A szakdolgozat beadható: ......................
1. fejezet Bevezetés Az Operációkutatás nevű tantárgy során rengeteg, korábban nem tárgyalt matematikai problémával és ezek megoldási módszereivel találkoztunk. Az egyik alapvető típusa ezeknek a függvények szélsőértékének kereséséhez kötődik. A vizsgált függvény jellegéből adódóan megkülönböztethetünk többek között lineáris, egész értékű és nemlineáris optimalizálást. Különösen az utóbbi terület az, ami kiemelkedik a többi közül - annak ellenére, hogy csak a múlt század közepén indult igazán fejlődésnek, mára nagyméretű problémák megoldását teszi lehetővé. A különböző programozási feladatok nagyon érdekesnek bizonyultak számomra. Azért választottam ezt a témát, mert úgy gondolom, hogy logikusan felépíthető és jól érthető a definíciók, tételek és módszerek sokasága. A szakdolgozatom két fő részből tevődik össze: az első az elmélet, a második a programozás. Így tehát egy rövid történeti áttekintés után ismertetem a feltételes és a feltétel nélküli optimalizálást, és az ezekhez kapcsolódó numerikus eljárásokat. Utóbbiak azért fontosak, mert a számítógépek számára kevésbé bonyolult, ezáltal gyorsabban értelmezhető algoritmusokat adnak, vagyis programozási szempontból előnyösebbek, mint a nem numerikus változatok. Ezután pedig a megvalósítás következik, mivel a választott témának része a web-es megjelenítés is, mely során az említett numerikus módszereket készítem el, böngészőben futtatható programként.
A kutató munka a Miskolci Egyetem stratégiai kutatási területén működő Fenntartható Természeti Erőforrás Gazdálkodás Kiválósági Központ / Alkalmazott Anyagtudomány és Nanotechnológia Kiválósági Központ / Mechatronikai és Logisztikai Kiválósági Központ / Innovációs Gépészeti Tervezés és Technológiák Kiválósági Központ keretében valósult meg.
5
2. fejezet Téma elméleti kifejtése 2.1. Történeti áttekintés A nemlineáris optimalizálás komoly múltra tekint vissza. Az első fontosabb mű, melyet meg kell említeni, Joseph-Louis Lagrange francia matematikus nevéhez fűződik. 1788ban készült el a Mécanique analitique (Analitikus mechanika) című könyve, melyben leírta a függvények szélsőértékeinek kiszámítására vonatkozó, egyenlőség feltételek esetén érvényes multiplikátoros módszerét. 1798-ban Lagrange kimondta a róla elnevezett egyenlőtlenségi elvet, amit később Osztrogradszkij és Gauss is megerősített. A következő említésre méltó állomás 1902, amikor Farkas Gyula magyar matematikus-fizikus publikálta a később róla elnevezett Farkas-lemmát (tételt). Ez a lineáris programozás egyik alaptétele, de a nemlineáris optimalizálás területén is használják. Az 1930-as években Gilbert Bliss és Constantin Carathéodory is egyenlőség feltételekkel rendelkező feladatokkal foglalkozott. Viszont 1939-ben William Karush az egyenlőtlenségi feltételekkel megadott problémákat vizsgálta, kutatásait azonban nem tette közzé. Utána nem sokkal Fritz John szintén ezt a témát kutatta, de eredményei gyengébbek voltak, mint Karushé. Az áttörést egy 1950-ben megjelent, Harold W. Kuhn és Albert W. Tucker által írt cikk hozta meg. Ők is azzal foglalkoztak, amivel Karush, és annak ellenére, hogy tőle függetlenül dolgoztak, mint utólag kiderült ugyanarra a megállapításra jutottak, mint Karush. Ebben az írásban szerepelt először a nemlineáris programozás, mint kifejezés. Sokan eddig nem ismerték fel a téma fontosságát és hasznosságát, de ekkor megindult a fejlődés ezen a téren. Egyre több helyen tudták felhasználni az újonnan elért eredményeket, például a mérnöki, valamint a közgazdasági tudományok területén.
2.2. Feltétel nélküli optimalizálás A kiinduló feladat a következő: adott egy f (x) : Rn → R függvény. Ennek kell a minimumát megkeresni. A feltétel nélküli jelző alapvetően azt jelenti, hogy nincs semmilyen korlátozás megadva az x döntési változóra.
6
2.2. Feltétel nélküli optimalizálás
2.2.1. Alapfogalmak Ahhoz, hogy a nemlineáris optimalizálási feladat megoldásához alkalmazandó módszert bemutassam, szükség van néhány definíció és tétel ismertetésére. Elsőként az iránymenti deriváltat, majd a differenciálhatóságot és a gradiens fogalmát definiáljuk. 2.1. definíció. Legyen S ⊆ Rn egy nemüres halmaz és legyen f : Rn → R függvény. Legyen x ∈ S és legyen d ∈ Rn , d ̸= 0 olyan vektor, hogy x + λd ∈ S, λ > 0 elegendően kicsi valós számra. Az f függvénynek az x pontban a d irány mentén vett iránymenti deriváltjának nevezzük és f ′ (x, d)-vel jelöljük az alábbi határértéket: f ′ (x, d) = lim
λ→0+
f (x + λ · d) − f (x) , λ
(amennyiben létezik a határérték.) Az iránymenti derivált a függvény d irányba történő változását mutatja. A differenciálhatóság a következőt jelenti: 2.2. definíció. Legyen S ⊆ Rn egy nemüres halmaz, és legyen f : Rn → R függvény. Az f függvényt differenciálhatónak mondjuk az S halmaz egy x belső pontjában, ha létezik egy ∇f (x) ⊆ Rn vektor és egy α : Rn → R függvény úgy, hogy f (x) = f (x) + ∇f (x)(x − x) + ∥x − x∥ α(x, x − x) minden x ∈ S-re, ahol limx→x α(x, x − x) = 0. Ennek van egy másik alakja, ahol x-et x-re, x-et pedig x + λd (λ ∈ R) vektorra cseréljük. Ebben az esetben x − x = λd. Így az előző összefüggés erre módosul: f (x + λd) = f (x) + λ∇f (x)d + |λ|∥d∥α(x, λd) bármely λ valós számra, és d vektorra. Ekkor a határérték: limλ→0 α(x, λd) = 0. A formulában szereplő ∇f (x) vektort hívjuk gradiensnek. Ehhez két tétel kapcsolódik. Az első az iránymenti derivált és a gradiens vektor közötti összefüggést mutatja meg. 2.3. tétel. Legyen f : S → R függvény differenciálható és legyen a d irányvektor egységhosszúságú, azaz |d| = 1. Ekkor f ′ (x, d) = ∇f (x)d. A második tétel pedig azt mutatja meg, hogy a gyakorlatban hogyan kell meghatározni egy függvény gradiensét. 2.4. tétel. Legyen f : S → R differenciálható függvény. A ∇f (x) vektor elemei az f függvény parciális deriváltjai, azaz ) ( ∂f (xn ) ∂f (x) ∂f (x) , ,..., . ∇f (x) = ∂x1 ∂x2 ∂xn A gradiens vektor előállításával analóg módon kapjuk a Hesse-mátrixot, melyre szintén szükség van a feladatok megoldásához. 7
2.2. Feltétel nélküli optimalizálás 2.5. definíció. Az f függvény második parciális deriváltjaiból álló H = ∇2 f (x) ∈ Rn×n mátrixot Hesse-mátrixnak nevezzük. (H(x))ij = (∇2 f (x))ij =
∂ 2 f (x) ∂xi ∂xj (
Mivel megegyeznek a másodrendű vegyes parciális deriváltak Hesse mátrix szimmetrikus.
∂ 2 f (x) ∂xj ∂xi
=
∂ 2 f (x) ∂xi ∂xj
) ,a
2.2.2. Az optimalitás szükséges és elégséges feltételei A minimum szükséges feltételére az alábbi tétel vonatkozik: 2.6. tétel. Legyen f : Rn → R függvény legalább egyszer differenciálható az x ∈ Rn pontban. Ekkor ha x lokális minimumpontja f -nek, akkor ∇f (x) = 0. Ehhez kötődik a stacionárius pont definíciója. 2.7. definíció. Az x ∈ Rn pontot stacionárius pontnak nevezzük, ha ∇f (x) = 0. Evidens, hogy a keresett minimumpont csak ezen pontok közül kerülhet ki. Amennyiben az f függvény legalább kétszer differenciálható, akkor a szükséges feltételt leíró 2.6 tétel így módosul: 2.8. tétel. Legyen f : Rn → R függvény legalább kétszer differenciálható az x ∈ Rn pontban. Ekkor ha x lokális minimumpontja f -nek, akkor ∇f (x) = 0 és H(x) Hesse mátrix pozitív szemidefinit. Ennek a tételnek a megfordítása egy kis módosítással pontosan az optimalitás elégséges feltételét adja meg. 2.9. tétel. Legyen f : Rn → R függvény legalább kétszer differenciálható az x ∈ Rn pontban. Ha ∇f (x) = 0 és H(x) pozitív definit, akkor x szigorúan lokális minimumpontja f -nek. Megjelent egy új tulajdonság, a definitség, amire a Hesse-mátrix vizsgálatakor van szükségünk. 2.10. definíció. Legyen A egy n × n-es (kvadratikus) mátrix. Az A mátrix pozitív definit, ha xT Ax > 0, ∀x ̸= 0, x ∈ Rn esetén. Ha megengedjük az egyenlőséget: xT Ax ≥ 0, akkor a mátrix pozitív szemidefinit. Analóg módon a mátrix negatív definit, ha xT Ax < 0, és negatív szemidefinit, ha xT Ax ≤ 0. Az előző kategóriákba nem sorolható mátrixokat indefinitnek nevezzük. Amennyiben f maximumának megtalálása a feladat, a fenti tételek ugyanígy érvényesek azzal a módosítással, hogy a H mátrixra vonatkozó részeknél pozitív helyett negatív definit a helyes. Egy konkrét feladat megoldásakor viszont célravezetőbb lehet úgy eljárni, hogy f maximuma helyett −f minimumát keressük.
8
2.2. Feltétel nélküli optimalizálás
2.2.3. A feladat megoldásának módja Mivel az f függvény globális minimumát szeretnénk megkapni, azt úgy oldhatjuk meg, hogy a lokális minimumpontokat megkeressük, majd kiválasztjuk közülük a legkisebb értékűt. Persze ilyenkor figyelni kell arra is, hogy a függvénynek van-e globális minimuma, tehát f korlátos legyen. A megoldás menete a szükséges és az elégséges tételek állításaiból következik. Csak olyan pontból lehet minimum, ami stacionárius, ezért először meg kell keresni az összes stacionárius pontot, tehát meg kell oldani a ∇f (x) = 0 egyenletrendszert. Ezután minden ilyen pontra megnézzük, hogy a Hesse mátrix pozitív definit-e az adott pontban. Ha igen, akkor a vizsgált pont lokális minimum. Ha van globális minimuma a függvénynek, akkor a lokális pontok közül a legkisebb értékű lesz a globális minimum pont. A feltétel nélküli optimalizálási feladatot tehát így oldjuk meg: 1. Kiszámítjuk azon x vektorok értékét, melyekre ∇f (x) = 0 teljesül, ezek eredményezik a lehetséges minimum helyeket. 2. Megvizsgáljuk a kapott x-ekre, hogy H(x) pozitív definit-e. Ha igen, akkor lokális minimumpontot kaptunk. 3. A függvény korlátosságát figyelembe véve, a lokális minimumpontok közül kiválaszthatjuk a globális pontot, ha van. Amennyiben a vizsgált mátrix szimmetrikus (tudjuk a Hesse mátrixról hogy az), vagyis HT = H, a pozitív definitség eldöntéséhez többféle mód áll rendelkezésünkre. • Gauss-elimináció használatával a vizsgálandó mátrixból egy felső háromszög mátrixot kapunk. Ha ennek főátlójában minden elem pozitív, akkor az eredeti mátrix pozitív definit. • A Hesse-mátrix minden k ×k-s főminor mátrixának kiszámolom a determinánsát. Ha ezek mindegyike pozitív, akkor a vizsgált mátrix pozitív definit. Ezt nevezik Főminor-módszernek. • Sajátérték teszt segítségével is megoldhatjuk a problémát - ha a H(x) minden sajátértéke pozitív, akkor a mátrix pozitív definit. (Mint tudjuk, a sajátérték az a λ szám, ami kielégíti az Hx = λx egyenletet) • Inercia tesztet használva a sajátértékek kiszámolása nélkül meg tudjuk állapítani, hogy a mátrixnak mennyi pozitív, mennyi zérus és mennyi negatív sajátértéke van. Ezt a három számot nevezzük inerciának. Ha csak pozitívak vannak, akkor a mátrix pozitív definit. A Cottle algoritmust kell használni ennél a módszernél. Példa Legyen f : R2 → R, és keressük f (x1 , x2 ) = x31 + x32 − 3x1 − 27x2 + 273 minimumát! ) ( 3x21 − 3 . Először is felírjuk a függvény gradiensét: ∇f (x) = 3x22 − 27 9
2.3. A feltétel nélküli optimalizálás numerikus módszerei Ezután meg } kell oldani a ∇f (x) = 0 egyenletrendszert.
3x21 − 3 = 0 3x22 − 27 = 0
. Gyors számolással megkapjuk, hogy x1 = ±1 és x2 = ±3. Összesen
tehát 4 olyan pontot kapunk, amelyek teljesítik a szükséges feltételt, vagyis stacionárius pontok: x1 = (1; 3), x2 = (1; −3), x3 = (−1; 3), x4 = (−1; −3). [ ] 6x1 0 A második lépés a Hesse mátrix felírása. H(x) = . Most pedig jöhet a 0 6x2 definitség[vizsgálata, a már meglévő pontokat a H mátrixba helyettesítve. ] 6 0 H(x1 ) = . Ez pozitív definit mátrix =⇒ x1 pont lokális minimumhely. 0 18 [ ] 6 0 H(x2 ) = . Ez indefinit mátrix =⇒ x2 pont ekkor nyeregpont. 0 −18 [ ] −6 0 H(x3 ) = . Ez is indefinit mátrix =⇒ x3 pont szintén nyeregpont. 0 18 ] [ −6 0 . Ez negatív definit mátrix =⇒ x4 pont lokális maximumhely. H(x4 ) = 0 −18 A függvény jellegéből adódóan nincs globális szélsőérték, csak lokálisak.
2.3. A feltétel nélküli optimalizálás numerikus módszerei Az előző alfejezetben láttuk, hogyan lehet a feltétel nélküli optimalizálás feladatát megoldani. Vannak azonban olyan módszerek erre a problémára, melyek a számítógépek számára jobban értelmezhetők, gyorsabb végrehajtást tesznek lehetővé. Ezek a numerikus metódusok, melyek nem feltétlenül fognak pontos eredményt szolgáltatni, ilyenkor viszont kellően kicsi hibahatár esetén nagyon közel lehetnek hozzájuk. Az alapfeladat nem változott, egy f (x) : Rn → R függvény minimumát kell meghatározni. Alapvetően 2 különböző típusú numerikus módszer létezik erre: a Newton jellegűek, valamint a vonalmenti optimalizálást alkalmazó eljárások. Szakdolgozatom témájából adódóan elsősorban az előbbi módszerekkel foglalkozom. Ezeken belül a továbbiakban a Newton módszert, valamint a kvázi Newton módszereket ismertetem. Később azonban a két típus együttes használatának lehetőségét is bemutatom.
2.3.1. Newton módszer A Newton módszer a numerikus analízisben többek között egyenletrendszerek eredményének kiszámítására alkalmazható. Azonban a nemlineáris optimalizálás esetén is lehet olyan formára hozni a feladatot, hogy ilyen módon legyen megoldható. Tehát először nézzük meg az egyenletrendszerek esetét. Legyenek a függvények f1 , f2 , . . . , fm : Rn → R többváltozósak és differenciálhatóak. Tegyük egyenlővé mindegyiket nullával! Az alábbi, m egyenletet tartalmazó rendszert kapjuk:
10
2.3. A feltétel nélküli optimalizálás numerikus módszerei
f1 (x) = 0 f2 (x) = 0 .. . fm (x) = 0 Alkalmazzuk az elsőfokú Taylor polinomot a függvények közelítésére egy xk ∈ Rn helyen. A következőt kapjuk: f1 (x) ≈ f1 (xk ) + ∇f1 (xk )(x − xk ) f2 (x) ≈ f2 (xk ) + ∇f2 (xk )(x − xk ) .. . fm (x) ≈ fm (xk ) + ∇fm (xk )(x − xk ) A mátrix-vektoros jelölésre áttérve, legyen F : Rn → Rm az a vektor, melynek elemei az f1 (x), f2 (x), . . . , fm (x) függvények. Emellett legyen egy J(xk ) mátrix, m × n dimenziós, melynek minden sorvektora az azonos sorszámú függvény xk pontbeli gradiense, vagyis J(xk ) = ∇F(xk ). Ez a J mátrix a Jacobi mátrix. A korábbi egyenletrendszer így átírható egyszerűbb alakra: F(x) ≈ F(xk ) + J(xk )(x − xk ). Az xk+1 -edik közelítést válasszuk meg úgy, hogy a Taylor polinom értéke nulla legyen, azaz F(xk ) + J(xk )(xk+1 − xk ) = 0, így xk+1 -et ez a lineáris egyenletrendszer adja meg: J(xk )(xk+1 − xk ) = −F(xk ). Ha az egyenletek száma és a változók száma megegyezik, akkor a J(xk ) mátrix négyzetes, és ha létezik inverze, akkor xk+1 = xk − J(xk )−1 F(xk ). A Newton módszer lényege tehát ez: egy tetszőleges x1 ∈ Rn vektorból kiindulva a J(xk )(xk+1 − xk ) = −F(xk ) , és xk+1 = xk − J(xk )−1 F(xk ) formulák segítségével meg lehet határozni az x2 , x3 , . . . , xk sorozatot, iteratív módon. Addig kell folytatni a számolást, amíg két, egymás utáni közelítés nem lesz elég közel, egy adott ε határon belül: ∥xk+1 − xk ∥ < ε. Nézzük, a feltétel nélküli optimalizálás esetén ez hogy működik. A cél megtalálni egy f : Rn → R, (egy n változós) függvény minimumát. Tudjuk, hogy a minimum szükséges feltétele a ∇f (x) = 0 , ami pontosan egy olyan egyenletrendszer, mely n 11
2.3. A feltétel nélküli optimalizálás numerikus módszerei változós és n egyenletből áll. Ennek megoldására éppen a Newton módszer szolgál. Ahhoz, hogy ez működjön, F(x) helyére ∇f (x)-et kell írni, így a Jacobi mátrix i-edik sora a ∇f (x) i-edik elemének a gradiense, azaz H(x) = ∇∇F(x) = ∇2 F(x). Ez a mátrix pedig éppen az f (x) függvény H(x) Hesse mátrixa. A Newton módszer így a következő: egy tetszőleges xk ∈ Rn (k = 1) kiindulási vektor használatával megoldjuk a H(xk )(xk+1 − xk ) = −∇f (xk ) egyenletrendszert, melyből meg tudjuk határozni a következő, xk+1 közelítést. Ha létezik H(x) inverze, a következő formula is használható: xk+1 = xk − H(xk )−1 ∇f (xk ). Célszerű bevezetni az sk = xk+1 − xk irányvektort. Az előző két formulát megkülönböztetjük elnevezésben, az újonnan bevezetett sk kiszámítására az alábbi képletek adódnak: Newton I. módszer: H(xk )sk = −∇f (xk ), Newton II. módszer:
sk = −H(xk )−1 ∇f (xk ).
Ha már ismerjük az sk irányvektort, a következő közelítést a xk+1 = xk + sk formulával kapjuk. Példa Newton módszerre Számítsuk ki Newton módszerrel az f (x1 , x2 ) = x21 + x1 x2 + 2x22 + x1 − 3x2 függvény minimumát az x1 = (2, 1) kiindulási vektort használva. Először is a gradiensvektort kell felírni, majd a Hesse mátrixot. ( ) [ ] 2x1 + x2 + 1 2 1 ∇f (x) = , H(x) = . x1 + 4x2 − 3 1 4 Ezekbe behelyettesítve az x1 pontot, azt kapjuk, hogy ( ) [ ] 6 2 1 ∇f (x1 ) = , H(x1 ) = . 3 1 4 A Newton I. módszert használva a H(xk )sk = −∇f (xk ) egyenletrendszert kell felírnunk, melyben s1 = (s1 , s2 ) irányvektorra a következő összefüggés adódik: 2s1 + s2 = −6 s1 + 4s2 = −3 12
2.3. A feltétel nélküli optimalizálás numerikus módszerei Ennek a megoldása hamar kiszámolható, az x1 -hez tartozó irányvektor ez lesz: s1 = (−3, 0). Így a következő(közelítés: x2 = x1 +s1 = (−1, 1). Ezután x2 behelyettesítésével foly) 0 tatjuk. ∇f (x2 ) = . Itt azonban megállunk, mert teljesül az optimalitás szükséges 0 feltétele, a gradiens értéke az adott pontban 0, tehát megkaptuk az optimumpontot, ami x2 = (−1, 1). A célfüggvény optimális értéke: fmin = −2. A Newton módszernél előfordulhatnak problémák, például a Hesse mátrix nem invertálható, vagy az irányvektor nem csökkenő irány (egyre távolodunk az optimumtól). Ezeket más módszerek használatával kerülhetjük el.
2.3.2. Kvázi Newton módszerek Ahogy a nevükből is sejteni lehet, a kvázi Newton módszerek a Newton módszer módosult változatai. A legfőbb eltérés az, hogy jelen esetben nem a Hesse mátrixot, vagy az inverzét használjuk, hanem ezeket közelítő mátrixokkal helyettesítjük. Legyen F : Rn → Rm függvény. Vizsgáljuk az F(x) = 0 egyenletrendszert. Ha az xk ∈ Rn pontbeli elsőfokú Taylor polinommal közelítjük a z F(x) függvényt, a függvényérték közelítése az xk+1 pontban F(xk+1 ) ≈ F(xk ) + J(xk )(xk+1 − xk ). Az optimalizálás esetén egy f : Rn → R függvény minimumának megtalálása a cél, alapvetően a ∇f (x) = 0 egyenletrendszert kell megoldani. Ezzel az előző összefüggés az alábbira módosul: ∇f (xk+1 ) ≈ ∇f (xk ) + H(xk )(xk+1 − xk ). Bevezetünk két új jelölést: yk a gradiensek különbsége, sk pedig az x vektorok különbsége, ekkor yk ≈ H(xk )sk . Most fogjuk a Hesse mátrixot helyettesíteni, mégpedig egy Bk közelítő mátrixszal. Ekkor yk ≈ Bk sk . A következő helyettesítést, a Bk+1 mátrixot úgy válasszuk meg, hogy yk = Bk+1 sk teljesüljön. Ezt az egyenlőséget nevezzük kvázi Newton feltételnek. Ha az összefüggést yk helyett sk -ra rendezzük, ezt kapjuk: sk ≈ H(xk )−1 yk . Ekkor, ha a Hesse mátrix inverzét helyettesítjük egy Dk közelítő mátrixszal: sk ≈ Dk yk . 13
2.3. A feltétel nélküli optimalizálás numerikus módszerei Dk+1 mátrixot pedig úgy kell megválasztani, hogy igaz legyen a következő: sk = Dk+1 yk Ezt szintén kvázi Newton feltételnek nevezzük. Mivel a Newton módszerből kiinduló eljárásokról van szó, alapvetően egy f : Rn → R függvény minimumát kell megtalálni, egy tetszőleges xk ∈ Rn vektorból kiindulva. Először a függvény gradiensét kell felírni. Ezután viszont a Hesse mátrixra már nincs szükség, mert helyettesítjük azt a korábban említett B-vel, vagy az inverzét a D-vel. A megfelelő egyenletrendszer megoldása után megkapjuk az sk vektort, majd ennek segítségével kiszámoljuk az xk+1 -et, és yk -t. xk+1 = xk + sk ,
yk = ∇f (xk+1 ) − ∇f (xk )
Végül a következő közelítő mátrixot is meg kell határozni. Ennek többféle módja lehetséges, ezek közül kettőt ismertetnék, melyek hasonló képleteket használnak. A BFGS eljárás A módszer a kidolgozóiról - Broyden, Fletcher, Goldfarb, Shanno - kapta nevét. A Hesse mátrixot egy szimmetrikus, pozitív definit B mátrixszal helyettesítjük. Az egyenletrendszer így a következő: Bk sk = −∇f (xk ), Az sk vektor, majd az xk+1 és az yk kiszámolása után ezzel a képlettel kell megadni Bk+1 -et, az újabb közelítő mátrixot: Bk+1 = Bk +
yk ykT (Bk sk )(Bk sk )T − . (Bk sk )T sk ykT sk
A DFP eljárás A módszer Davidon, Fletcher, és Powell kutatókról kapta nevét. Itt a Hesse mátrix inverzét helyettesítjük egy szimmetrikus, pozitív definit D mátrixszal. Az egyenletrendszer a következő lesz: sk = −Dk ∇f (xk ), Ebből kiszámoljuk az sk vektort, majd az xk+1 és az yk vektorokat. Utána az újabb közelítő mátrixot, Dk+1 -et ezzel a képlettel kell meghatározni: Dk+1 = Dk +
(Dk yk )(Dk yk )T sk sT k − . (Dk yk )T yk sT k yk
Az eljárások kilépési feltételeként használhatunk egy ε pontossági értéket - ha ∥xk+1 − xk ∥ < ε, akkor megállunk. 14
2.3. A feltétel nélküli optimalizálás numerikus módszerei Példa kvázi Newton módszerre Oldjuk meg a Newton módszernél már korábban szereplő feladatot, DFP eljárással. f (x1 , x2 ) = x21 + x1 x2 + 2x22 + x1 − 3x2
függvény minimumát keressük, [ ] 1 1 az x1 = (2, 1) kiindulási vektort használva. Legyen a D1 = 1 2 Először felírjuk a gradiensvektort, majd kiszámoljuk értékét az x1 pontban. ( ) 2x1 + x2 + 1 ∇f (x) = , ∇f (x1 ) = (6, 3). x1 + 4x2 − 3 Ezután meghatározzuk s1 értékét, majd ezt felhasználva x2 -t. [ ][ ] [ ] 1 1 6 −9 s1 = −D1 ∇f (x1 ) = − = 1 2 3 −12 ) ) ( ( ) ( −7 −9 2 = + x2 = x1 + s1 = −11 −12 1 Most jön a D2 mátrix kiszámítása. Ehhez a következőkre lesz szükség: ∇f (x2 ) = (−24, −54) y1 = ∇f (x2 ) − ∇f (x1 ) = (−30, −57) ] ] [ ] [ [ [ ] −30 ] 81 108 −9 [ T T = 954 s1 y1 = −9 −12 s1 s1 = −9 −12 = −57 108 144 −12 ] ] [ ] [ ][ [ [ ] −30 −87 1 1 −30 = 10818 (D1 y1 )T y1 = −87 −144 = D1 y1 = −57 −144 1 2 −57 [ ] ( ) [ ] 7569 12528 −87 (D1 y1 )(D1 y1 )T = −87 −144 = 12528 20736 −144 Ezeket felhasználva a D2 mátrix: [ ] [ ] [ ] [ ] 1 1 81 108 7569 12528 0, 3852 −0, 0449 1 1 D2 = + − = 954 108 144 10818 12528 20736 1 2 −0, 0449 0, 2341 Így épül fel a DFP kvázi Newton módszer, a folytatásban ki kéne számolni az s2 = −D2 ∇f (x2 ) vektort, innen ismételni a fenti lépéseket, amíg a kilépési feltétel nem teljesül egy adott ε-ra. A fejezet elején említettem, hogy a Newton jellegű módszerek és a vonalmenti keresést alkalmazó eljárás együttesen is használhatóak, ennek bemutatásához az utóbbi típust röviden ismertetem, majd leírom az egyesített algoritmusok működési elvét.
15
2.3. A feltétel nélküli optimalizálás numerikus módszerei
2.3.3. Newton módszerek vonalmenti minimumkereséssel A vonalmenti keresés lényege a következő: egy x1 ∈ Rn pontból indulunk ki, majd választunk egy d1 ∈ Rn irányvektort. A d irányban megkeressük a függvény minimumát. Ez egy egyváltozós minimalizálási feladat. A φ(λ) = f (xk + λdk ) tehát egyváltozós függvény, ennek legyen a minimuma λk . Ezt felhasználva az x vektor következő közelítése: xk+1 = xk + λk dk A vonalmenti keresést alkalmazó módszerek alapvetően abban különböznek egymástól, hogy a kereső irányt hogyan határozzák meg. Ettől függetlenül viszont mindegyiknél szükség van egyváltozós függvény minimumának kiszámítására. Egyváltozós optimalizálás Az egyváltozós minimalizálási feladat: f (x) → min, ahol x ∈ R. Minden numerikus eljárás, ami ezt megoldja, alapvetően bizonytalansági intervallumot használ. Ez egy olyan [a, b] intervallum, amit mi határozunk meg, és tartalmazza a minimumpontot. A tartományt az algoritmusok folyamatosan szűkítik, mégpedig osztópontok segítségével. Ezekben az osztópontokban (általában kettő darab) ki kell számolni a függvény értékét, majd ez alapján meg lehet határozni az új bizonytalansági intervallumot. Ennek módját az alábbi tétel írja le: 2.11. tétel. Legyen f : R → R szigorúan kvázikonvex függvény az [a, b] bizonytalansági intervallumon, továbbá c, d ∈ (a, b) úgy, hogy c < d. Ekkor Ha f (c) > f (d) ⇒ f (x) = f (d) minden x ∈ [a, c) esetén ⇒ új intervallum: [c, b]. Ha f (c) 5 f (d) ⇒ f (x) = f (c) minden x ∈ (d, b] esetén ⇒ új intervallum: [a, d]. A függvény konvexitási tulajdonsága szükséges ahhoz, hogy az eljárás konvergens legyen. Ha nem vizsgáljuk meg a konvexitást, akkor nincs garancia az eljárások konvergenciájára. Az eltérő módszerek leginkább abban különböznek, hogy a belső pontokat hogyan veszik fel, számolják ki. Számos módszer létezik, például az uniform keresés, a Dichotomous keresés, a Fibonacci keresés, a Bisection módszer. Én a későbbi feladatmegoldáshoz a Golden section, vagyis aranymetszés keresést fogom alkalmazni, mivel az a leggyakrabban használt eljárás. A Golden section egy olyan módszer, melynél az új intervallumok mindkét esetben ugyanolyan hosszúak lesznek. Az osztópontokat úgy kell megválasztani, hogy az aktuális [a, b] hosszának (ezt jelölje L,) α-szorosai legyenek az új tartományok, tehát b − c = αL, és d − a = αL. Ekkor a következő képletek adják meg az osztópontokat: ck = bk − αLk = ak + (1 − α)Lk , dk = ak + αLk . α értéke szabadon választható az ( 12 , 1) tartományban, de célszerű úgy választani, hogy minél kevesebbszer kelljen a függvény értékét kiszámolni. Ilyen lépést úgy lehet megspórolni, hogyha az új bizonytalansági intervallum egyik belső pontja megegyezik az előző intervallum egyik belső pontjával. Ekkor ugyanis az egyik pontban lévő
16
2.3. A feltétel nélküli optimalizálás numerikus módszerei függvényértéket már egyszer meghatároztuk, nem kell újra. Ez azt jelenti, hogy az [ak+1 , bk+1 ] = [ck , bk ] intervallumban: ck+1 ak + (1 − α)Lk + (1 − α)αLk ak + (1 − α)Lk (1 + α) (1 − α)Lk (1 + α) (1 − α)(1 + α)
= = = = =
dk , vagyis ak + αLk , ennek levezetése: ak + αLk αLk α
A fenti egyenlőség tehát pontosan akkor teljesül, ha α2 + α − 1 = 0. Az [ak+1 , bk+1 ] = [ak , dk ] bizonytalansági intervallum esetén analóg módon szintén ezt az egyenletet kapjuk. Ennek megoldása az aranymetszésnél, mint fogalomnál általáno√ 5−1 san ismert α = 2 ≈ 0, 618 érték. A Golden section eljárás során érdemes ezt az α értéket használni. Összefoglalva tehát a módszert: Ha f (ck ) > f (dk ) ⇒ [ak+1 , bk+1 ] = [ck , bk ], ck+1 = dk , dk+1 = ak+1 + αLk+1 Ha f (ck ) 5 f (dk ) ⇒ [ak+1 , bk+1 ] = [ak , dk ], ck+1 = ak+1 + (1 − α)Lk+1 , dk+1 = ck . Az eljárást egészen addig csináljuk, amíg az intervallum hossza nagyobb, mint egy ε hibahatár kétszerese. Ha ugyanis kisebb lesz ennél a hossz, akkor az intervallum felezőpontját véve minimumnak, a hibahatáron belül lesz ez a közelítés. Tehát a kilépési feltétel: Lk < 2ε, valamilyen k-ra. Most, hogy az egyváltozós optimalizálás egy metódusát bemutattam, visszatérünk a többváltozós függvényekhez, és leírom azt, hogy erre egyáltalán miért volt szükség. Tehát a Newton típusú módszereknél úgy kapjuk meg a következő közelítést az x vektorra, hogy egyszerűen hozzáadjuk a számolás során kapott s vektort. Ezen a ponton azonban lehet vonalmenti minimumkeresést is végezni. Ilyenkor az s-t iránynak tekintjük, és ebben az irányban elvégzünk egy vonalmenti optimalizálást. A lépések tehát a következők: az irányvektort az s adja meg, dk = sk . A dk irányban megkeressük a függvény minimumát, vagyis megoldjuk a φ(λ) = f (xk + λdk ) → min! egyváltozós minimalizálási feladatot (például a Golden section eljárással), legyen a kapott érték λk . Ezt felhasználva az x vektor következő közelítése: xk+1 = xk + λk dk Ezt a vonalmenti keresést természetesen minden Newton típusú eljárással, így a Newton módszer mellett a kvázi Newton módszerekkel is lehet alkalmazni. Ezzel a feltétel nélküli optimalizálás témaköre véget ért, ismertettem az alapokat, a feladat megoldásának módját, és a kapcsolódó numerikus módszerek közül bemutattam a Newton jellegűeket, a továbbiakban a feltételes optimalizálásról lesz szó.
17
2.4. Feltételes optimalizálás
2.4. Feltételes optimalizálás 2.4.1. Alapfogalmak Legyen Df az f (x) függvény értelmezési tartománya, és X egy nyílt halmaz. Egy optimalizálási feladatot akkor nevezünk feltételesnek, ha az x ∈ Rn döntési változóra az x ∈ Df ⊂ Rn és az x ∈ X ⊂ Rn előírások mellett további kikötéseket teszünk. S-sel jelöljük a feltételeket meghatározó halmazt, és megengedett (lehetséges) megoldások halmazának vagy feltételi halmaznak nevezzük. A feladat tehát a következő: meg kell találnunk az S halmaz azon pontját, amelyben az f : Rn → R célfüggvény értéke minimális. A feltételes optimalizálási feladat jelölése: min{f (x) : x ∈ S} A feladat megoldásához az alábbi fogalmakat használjuk kiindulásként: 2.12. definíció. Legyen az f : Rn → R függvény. Az f függvény x ∈ Rn pontbeli javító irányok kúpjának nevezzük és F -fel jelöljük az alábbiakban definiált halmazt: F = {d : d ∈ Rn , d ̸= 0, f (x + λd) < f (x) minden λ ∈ (0; δ) -ra δ > 0 esetén}. Ez azt jelenti, hogy ha x pontból legalább egy kicsit, folyamatosan elmozdulva a d irányba, a célfüggvény értéke javul, akkor a d javító irány, így d ∈ F . 2.13. definíció. Legyen az S ⊂ Rn nemüres halmaz és x az S halmaz lezártjának ey pontja. Az S halmaz x pontbeli lehetséges irányok kúpjának nevezzük és D-fel jelöljük az alábbiakban definiált halmazt: D = {d : d ∈ Rn , d ̸= 0, x + λd ∈ S minden λ ∈ (0; δ) -ra δ > 0 esetén} Ez arról szól, hogy ha x pontból legalább egy kicsit, folyamatosan elmozdulva a d irányba, benne maradunk az S tartományban, akkor d lehetséges irány, vagyis d ∈ D. Ennek a két halmaznak a kapcsolatát mutatja meg a következő tétel, a feltételes optimalizálás geometriailag szükséges feltétele. 2.14. tétel. Legyen az S ⊂ Rn nemüres halmaz és f : S → R függvény. Adott a min{f (x) : x ∈ S} feltételes optimalizálási feladat. Ha x ∈ S lokális minimumpont, akkor F ∩ D = ∅. A tétel értelmezése: lokális minimumpont esetén nem létezik olyan irány, amerre ha elmozdulnánk, akkor javulna a célfüggvény és a tartományban is benne maradnánk. Bizonyítás. Indirekt módon tegyük fel, hogy F ∩ D ̸= ∅. Ekkor van olyan d ̸= 0 vektor, amelyre d ∈ F és d ∈ D. Tehát a d vektor javító irány és lehetséges irány is egyben, ami ellentmond annak, hogy az x pont optimális megoldás. A továbbiakban feltesszük, hogy a célfüggvény differenciálható, valamint az S halmazt függvények segítségével adjuk meg. Ezekkel a függvényekkel szembeni előírás a folytonosság és a differenciálhatóság. Így megfogalmazhatjuk a javító irányok kúpjának és a lehetséges irányok kúpjának karakterizációját.
18
2.4. Feltételes optimalizálás Legyen az S feltételi halmaz a következő: S = {x : x ∈ X, gi (x) 5 0, i = 1, 2, . . . , m}, ahol g1 , g2 , . . . , gm : Rn → R függvények és X ⊆ Rn nemüres, nyílt halmaz. A gi (x) 5 0 jelzi, hogy egyenlőtlenségi feltételekről van szó. Ha egyenlőségi előírások is vannak, akkor a tételek egy kicsit változnak, bővülnek, de a módszer bemutatása szempontjából elegendő a csak egyenlőtlenségi feltételek esetével foglalkozni. Fontos még megjegyezni, hogy ha a célfüggvény minimumát keressük, és az egyenlőtlenségi feltételeknél 5 jel szerepel, akkor standard optimalizálási feladatról beszélünk, a továbbiakban csak erről lesz szó. (A nem standard feladatokat köztes lépésekkel, átalakításokkal általában visszavezetik standard alakra.) 2.15. definíció. Tekintsük az S = {x : x ∈ X, gi (x) 5 0, i = 1, 2, . . . , m} feltételi halmazt, és egy x ∈ S pontot. A gi (x) függvény aktív a x pontban, ha gi (x) = 0 és inaktív, ha gi (x) < 0. 2.16. tétel. Legyen f : Rn → R függvény differenciálható x ∈ Rn pontban. Ekkor az F0 = {d : ∇f (x)d < 0} halmaz javító irányok kúpja, és F0 ⊆ F . 2.17. tétel. Legyen S = {x : x ∈ X, gi (x) 5 0, i = 1, 2, . . . , m}, ahol g1 , g2 , . . . , gm : Rn → R függvények és X ⊆ Rn nemüres, nyílt halmaz, továbbá adott egy x ∈ S pont. Legyen I = {i : gi (x) = 0} az aktív feltételi függvények indexhalmaza és tegyük fel, hogy a gi , i ∈ I függvények az x pontban differenciálhatók és a gi , i ̸∈ I függvények az x pontban folytonosak. Ekkor a G0 = {d : ∇gi (x)d < 0, i ∈ I} halmaz az S halmaz x pontbeli lehetséges irányok kúpja, és G0 ⊆ D. Az előző két tételhez kapcsolódik a geometriai-algebrai szükséges feltétel, mely szerint (a tételek feltevései és jelölései mellett): Ha x lokális minimumpont, akkor F0 ∩ G0 = ∅. Szintén a feltételes optimalizálás szükséges feltételét írja le a Fritz John tétel: 2.18. tétel. Tekintsük a min{f (x) : x ∈ S} feltételes optimalizálási feladatot, ahol S = {x : x ∈ X, gi (x) 5 0, i = 1, 2, . . . , m}, amelynél f, g1 , g2 , . . . , gm : Rn → R függvények és X ⊆ Rn nemüres, nyílt halmaz. Legyen adott továbbá egy x ∈ S pont és legyen I = {i : gi (x) = 0} indexhalmaz és tegyük fel, hogy az f, gi , i ∈ I függvények az x pontban differenciálhatók, valamint a gi , i ̸∈ I függvények az x pontban folytonosak. Ha x lokális optimális megoldás, akkor léteznek u0 , ui , i ∈ I számok úgy, hogy ∑ ui ∇gi (x) = 0 u0 ∇f (x) + i∈I
(u0 , uI ) (u0 , uI )
= ̸=
0 0 19
2.4. Feltételes optimalizálás ahol uI egy olyan vektor, melynek komponensei ui = 0, i ∈ I. De feltehetjük azt is, hogy f, gi , i = 1, 2, . . . , m függvények az x pontban differenciálhatók (tehát nem csak i ∈ I esetén), akkor a tétel így módosul: Ha x lokális optimális megoldás, akkor léteznek u0 , ui , i ∈ I számok úgy, hogy u0 ∇f (x) +
m ∑
ui ∇gi (x)
=
0
ui gi (x) (u0 , u) (u0 , u)
= = ̸=
0, 0 0
i=1
i = 1, 2, . . . , m
ahol u egy olyan vektor, melynek komponensei ui = 0, i = 1, 2, . . . , m. A gyakorlatban ez azt jelenti, hogy meg lehet adni olyan, Lagrange szorzóknak nevezett u0 , ui számokat, hogy a függvények gradienseinek az u0 , ui számokkal vett lineáris kombinációi eredményként nullát adjanak. A tétel második verziójába bekerült egy új feltétel is. Ez azért van, mert az elsőnél csak az aktív feltételeket nézzük, amikre definíció alapján teljesül a gi (x) = 0, tehát nem kell kikötni, hogy ui gi (x) = 0 legyen. Viszont az inaktív feltételekre ez nem igaz, tehát az utóbbi összefüggés is megjelenik a tételben. Ez egyben azt is eredményezi, hogy az összes inaktív feltételnél csak akkor teljesül ez az egyenlőség, ha a hozzájuk tartozó ui = 0.
2.4.2. KKT feladat Ezeket a Fritz John optimalitási feltételeket használva nem egyértelmű a Lagrange szorzók meghatározása. A továbbiakban ezért olyan Lagrange szorzókkal dolgozunk, melyeknél a célfüggvény együtthatója 1. Ezt az új feltételrendszert a Karush-KuhnTucker tétel írja le, erre épül az optimalizálási feladat megoldása. 2.19. tétel. Tekintsük a min{f (x) : x ∈ S} feltételes optimalizálási feladatot, ahol S = {x : x ∈ X, gi (x) 5 0, i = 1, 2, . . . , m}, amelynél f, g1 , g2 , . . . , gm : Rn → R függvények és X ⊆ Rn nemüres, nyílt halmaz. Legyen adott egy x ∈ S pont és legyen I = {i : gi (x) = 0} indexhalmaz és tegyük fel, hogy az f, gi , i ∈ I függvények az x pontban differenciálhatók, valamint a gi , i ̸∈ I függvények az x pontban folytonosak. Tegyük fel továbbá, hogy a ∇gi (x), i ∈ I vektorok lineárisan függetlenek. Ha x lokális optimális megoldás, akkor léteznek egyértelműen ui , i ∈ I számok úgy, hogy ∑ ∇f (x) + ui ∇gi (x) = 0 i∈I
ui
=
0,
i∈I
Viszont feltehetjük azt is, hogy f, gi , i = 1, 2, . . . , m függvények az x pontban differenciálhatók (tehát nem csak i ∈ I esetén), akkor ez lesz a tétel:
20
2.4. Feltételes optimalizálás Ha x lokális opt. megoldás, akkor léteznek egyértelműen ui , i ∈ I számok úgy, hogy ∇f (x) +
m ∑
ui ∇gi (x)
=
0
ui gi (x) ui
= =
0, 0,
i=1
i = 1, 2, . . . , m i = 1, 2, . . . , m
A Karush-Kuhn-Tucker tétel alapján három feltételcsoportot különböztetünk meg egymástól: • PF: primál feltételek Az x ∈ S lehetséges megoldásra vonatkozó gi (x) 5 0, i = 1, 2, . . . , m; x ∈ X feltételek. • DF: duál feltételek A tételben szereplő összefüggések: ∇f (x) +
m ∑
ui ∇gi (x)
=
0
ui
=
0,
i=1
i = 1, 2, . . . , m
• KF: komplementaritási feltételek A tételben lévő ui gi (x) = 0 i = 1, 2, . . . , m egyenlőségek. A feltételeket együttesen Karush-Kuhn-Tucker (KKT) optimalitási feltételeknek nevezzük. Ha egy x ponthoz léteznek olyan u = (u1 , u2 , . . . , um ) Lagrange szorzók, melyekre az (x, u) vektor a KKT feltételeket kielégíti, az x pontot KKT pontnak hívjuk. Ilyenkor a tételben szereplő ∇gi (x), i ∈ I feltevés miatt a Lagrange szorzók egyértelműen meghatározhatók. Egy feladat optimális megoldása csakis KKT pont lehet, mert a szükséges feltételeknek azok felelnek meg, viszont ez még nem garantálja, hogy az adott pont az optimum, ráadásul több KKT pontot is kaphatunk egy feladaton belül. Éppen ezért az optimalitás szükséges feltételei mellett a feladatok megoldásához az elégséges feltételek is kellenek. Ezek bemutatásához először a Lagrange függvényt definiálom, majd annak aktualizált változatát, valamint a reguláris pontot. A meghatározások az egyenlőséges feltételeket is tartalmazó, tehát általánosabb optimalizálási feladatra vonatkoznak, de a v Lagrange szorzókat és a hi (x) egyenlőséges feltételeket elhagyva a csak egyenlőtlenséges feladatnak megfelelő képleteket kapnánk. 2.20. definíció. A Lagrange függvény az optimalizálási feladatban szereplő célfüggvény, az egyenlőtlenséges és az egyenlőséges feltételi függvények lineáris kombinációja, amelyben a célfüggvény együtthatója 1. A feltételi függvények együtthatói a már korábban említett Lagrange szorzók. A Lagrange függvény képlete: L(x, u, v) = f (x) +
m ∑ i=1
ui gi (x) +
k ∑
vi hi (x)
i=1
21
2.4. Feltételes optimalizálás 2.21. definíció. Legyen x vektor az optimalizálási feladat KKT pontja és legyenek u, v vektorok a KKT feltételeket kielégítő Lagrange szorzók. Ezekhez az u, v Lagrange szorzókhoz tartozó Lagrange függvényt aktualizált Lagrange függvénynek nevezzük, melynek képlete: La (x, u, v) = f (x) +
m ∑
ui gi (x) +
i=1
k ∑
v i hi (x)
i=1
Az La függvényben minden egyenlőtlenséges feltétel szerepel, később a feladatmegoldás során ki kell szűrni az aktívakat. 2.22. definíció. Az x ∈ S lehetséges megoldást regulárisnak nevezzük, ha az aktív feltételekhez tartozó ∇gi (x), i ∈ I vektorok lineárisan függetlenek. Ezekkel a definíciókkal már meg lehet fogalmazni a feltételes optimalizálási feladat elégséges feltételét. Van gyenge és erős feltétel, valamint a leggyakrabban használható, én utóbbit ismertetem. 2.23. tétel. Legyen X ⊆ Rn nemüres nyílt halmaz, f, g1 , g2 , . . . , gm , h1 , h2 , . . . , hk : Rn → R kétszer differenciálható függvények. Tekintsük az f (x) → min! gi (x) 5 0
i = 1, 2, . . . , m
hi (x) = 0 i = 1, 2, . . . , k x∈X matematikai programozási feladatot. Jelölje S a feltételi halmazt. Legyen x vektor az optimalizálási feladat reguláris KKT pontja és legyenek u, v vektorok a KKT feltételeket kielégítő Lagrange szorzók. Legyen az aktualizált Lagrange függvény Hesse mátrixa az x = x pontban. Az I = {i : gi (x) = 0} indexhalmaz legyen felosztva az I + = {i ∈ I : ui > 0} és az I 0 = {i ∈ I : ui = 0} indexhalmazokra. Definiáljuk az alábbi A alteret A = {d ̸= 0 : ∇gi (x)d = 0 i ∈ I + , ∇hi (x)d = 0 i = 1, 2, . . . , k}. Ha dT ∇2 La (x)d > 0 minden d ∈ A irányra, akkor az x KKT pont szigorú lokális minimumpontja az optimalizálási feladatnak. Egy feladat megoldása során a fenti feltételt a következő módon használjuk: ha már megvan egy KKT pont, akkor felírjuk a Lagrange függvényt, majd az adott KKT pontra vonatkozó aktualizált Lagrange függvényt is. Ezután létrehozunk egy C mátrixot, ilyen felépítéssel: [ ] [ ] + ∇g (x), i ∈ I , (u > 0) A BT i i , C= , ahol A = ∇2 La (x), és B = ∇hi (x), i = 1, 2, . . . , k B 0 tehát az A mátrix az aktualizált Lagrange függvény Hesse mátrixa, míg a B mátrix sorvektorai az aktív egyenlőtlenséges feltételi függvények és az egyenlőséges feltételi függvények x pontbeli gradiensei. Ha ez a C mátrix pozitív definit, akkor az x pont az optimalizálási feladat szigorú lokális minimumpontja. 22
2.4. Feltételes optimalizálás Példa feltételes optimalizálási feladatra Oldjuk meg az alábbi, egyenlőtlenséges feltételes optimalizálási feladatot! f (x) = x21 + x22 → min! g1 (x) = 2 − x1 − x2 5 0 g2 (x) = −x1 + x2 − 1 5 0 X = R2 Első lépésként ( a)függvények gradienseit. ( ) ( ) meghatározzuk −1 −1 2x1 , ∇g1 (x) = , ∇g2 (x) = . ∇f (x) = −1 1 2x2 ( ) ( ) ( ) −1 −1 2x1 + u2 . A gradiens vektorokra vonatkozó lineáris kombináció: + u1 1 −1 2x2 Ezek segítségével felírjuk a KKT optimalitási feltételeket: 2x1 − u1 − u2 = 0 2x2 − u1 + u2 = 0 DF u ,u = 0 1
2
}
u1 (2 − x1 − x2 ) = 0
KF
u2 (−x1 + x2 − 1) = 0 } 2 − x1 − x2 5 0
PF
−x1 + x2 − 1 5 0
Meg kell nézni, vannak-e KKT pontok. Az u1 , u2 = 0 feltétel miatt 4 esetet kell megvizsgálnunk: u1 = 0, u2 = 0; u1 > 0, u2 = 0; u1 = 0, u2 > 0 vagy u1 > 0, u2 > 0. A KKT feltételekkel összevetve a fenti kikötéseket, egyedül a második esetben adódik KKT pont. Itt x1 = x2 = 1, és u2 = 0 mellett u1 = 2. Ezután meg kell vizsgálni, hogy a talált KKT pont optimális megoldás-e. A Lagrange függvény: L(x1 , x2 , u1 , u2 ) = x21 + x22 + u1 (2 − x1 − x2 ) + u2 (−x1 + x2 − 1) 2x1 − u1 − u2 2x2 − u1 + u2 ∇L = 2 − x1 − x2 −x1 + x2 − 1
2
0
−1
1
−1 −1
0 2 −1 ∇ L= −1 −1 0 2
0
1 0 0
Aktualizáljuk a Lagrange függvényt az u vektorral(- u1 = 2,)u2 = 0. [ ] 2x − 2 2 0 1 ∇La = ∇2 La = La (x1 , x2 ) = x21 + x22 + 2 (2 − x1 − x2 ) 2x2 − 2 0 2
23
2.5. A feltételes optimalizálás numerikus módszerei
Az x = (1; 1) pontban ∇2 La (x) =
[ ] 2 0
= az A mátrix. A B mátrixba csak ∇g1 [ ] kerül, mert u2 = 0 miatt a g2 inaktív feltétel. Így B = −1 −1 . A C mátrix pedig: 2 0 −1 Erről (pl. inercia teszttel) belátható, hogy pozitív definit, ezért C= 0 2 −1 −1 −1 0 az x = (1; 1) KKT pont szigorú lokális minimumpont. Ezzel a feladat kész. 0 2
2.5. A feltételes optimalizálás numerikus módszerei Az előző alfejezetben láthattuk, hogyan lehet megoldani a feltételes optimalizálási feladatot. A KKT módszer azonban programozási szempontból nem ideális, ezért több olyan, közelítést alkalmazó numerikus eljárás létezik, ami szintén alkalmasak a feladat elvégzésére. Alapvetően két különböző módszercsalád van. Az első a lehetséges irányok módszere. Ez a feltétel nélküli optimalizálás során használható vonalmenti minimalizálásra épül, és lépésenként keresi a feladat célfüggvényének optimális megoldását egy-egy irányban. Ilyen eljárás a Zoutendijk módszer, illetve a Rosen-féle gradiens projekciós módszer. A másik típusú módszercsalád lényege az, hogy a feltételes opt. feladatot visszavezetjük feltétel nélküli feladatok ismételt megoldására. Ez az eljárás a SUMT (Sequential Unconstrained Minimization Technique) nevet kapta. Működésének lényege: beépíti a feltételeket a célfüggvénybe, így egy olyan új feladatot kapunk, mely már feltétel nélküli, és ezt oldjuk meg sorozatosan más-más paraméterrel. Két változata van: a korlátfüggvényes és a büntetőfüggvényes eljárás. Előbbi esetén ún. korlátfüggvény segítségével lehet beépíteni a célfüggvénybe a feltételeket, ezzel elérve, hogy a lépések során soha nem lépünk ki a megengedett megoldások tartományából. Utóbbi fajtát fogom részletesen bemutatni a következőkben.
2.5.1. Büntetőfüggvényes SUMT A SUMT módszer tehát feltétel nélküli feladatokra vezeti vissza a feltételes opt. feladatot. Tulajdonképpen az történik, hogy a feltételi függvényeket egy büntetőfüggvényre cserélve hozzáadjuk a célfüggvényhez, majd az új célfüggvényt - nevezzük segédfüggvénynek - kell feltétel nélküli módszerrel megoldani. A büntetőfüggvénynek van egy paramétere, melyet növelve a segédfüggvény minimuma egyre közelebb kerül az eredeti célfüggvény minimumához. Ezeket a jelöléseket használjuk: a segédfüggvény - ϕ tehát az új célfüggvényt jelenti. A büntető paraméter egy µ > 0 szám, amely a büntetőfüggvénnyel (B) együtt alkotja a büntetőtagot (µB). Az új, feltétel nélküli opt. feladat alakja ez lesz: ϕ(x) = f (x) + µB(gi (x), hi (x)) → min! x ∈ Rn 24
2.5. A feltételes optimalizálás numerikus módszerei A büntetőfüggvénynek olyan folytonos függvényt kell választani, ami pozitív büntető értéket ad, ha nem teljesednek a feltételek, viszont teljesülés esetén zérus legyen a büntetés. Természetesen több ilyen függvényt lehet találni, de legtöbbször egyenlőséges feltételnél az y 2 függvényt, és egyenlőtlenséges esetben a [max{0, y}]2 függvényt választjuk. Leggyakrabban tehát egy általános alakú feltételes opt. feladatnál f (x) → min! gi (x) 5 0 i = 1, 2, . . . , m hi (x) = 0 i = 1, 2, . . . , k x∈X az alábbi segédfüggvény kerül használatba: [ m ] k ∑ ∑ [max{0, gi (x)}]2 + [hi (x)]2 . ϕ(x) = f (x) + µ i=1
i=1
A módszer során tehát ezt a ϕ(x) → min! feltétel nélküli opt. feladatot oldjuk meg lépésenként növelve a µ büntető paraméter értékét. Egy tetszőleges x ∈ Rn vektorral, valamint általános gyakorlat szerint µ = 1 értékkel kezdjük el. Később minden új lépés során µ-t (ez szintén szokás) 10-szeresére növeljük. Ha a ϕ(x) = min{f (x) + µk B(x) : x ∈ X} minimalizálási feladat optimális megoldása xk+1 , akkor megállási feltétel lehet a µk B(xk+1 ) < ε teljesülése.
Példa büntetőfüggvényes SUMT eljárásra Oldjuk meg a következő optimalizálási feladatot: f (x) = (x1 − 1)2 + x22
→ min!
h(x) = x1 − x2 = 0 Itt nincs egyenlőtlenségi feltétel, egyenlőségi pedig egy van. Így a feltétel nélküli optimalizálási segédfeladat az alábbi: ϕ(x) = (x1 − 1)2 + x22 + µ(x1 − x2 )2
→ min!
Legyen a kiinduló vektor az x = (5, 10), a µ büntető paraméter pedig az általános szokás szerint 1, minden lépésnél 10-szeresére növelve. Az eredmények a táblázatban láthatóak. µ
xµ
ϕ(xµ )
f (xµ )
1
(0.666667; 0.333333)
0.33333333
0.22222222
10
(0.52381; 0.47619)
0.47619048
0.45351474
100
(0.502488; 0.497512)
0.49751244
0.49503725
1000
(0.50025; 0.49975)
0.49975012
0.49950037
10000
(0.500025; 0.499975)
0.499975
0.49995
100000
(0.500002; 0.499998)
0.4999975
0.499995 25
2.5. A feltételes optimalizálás numerikus módszerei Az adatokból leolvasható, hogy a korábban említett összefüggéseknek megfelelően a µ növelésével egyre közelebb kerülünk az optimumhoz, a segédfüggvény és az eredeti f függvény értékei is a kapott helyen egyre kisebb eltérést mutatnak. Belátható, hogy a megoldás: x = (0.5, 0.5), és a célfüggvény értéke ebben a pontban szintén 0.5 lesz. Ennek a módszernek a bemutatásával zárult az elméleti összefoglalás, a következő fejezet a program megvalósítását és annak részleteit taglalja.
26
3. fejezet Fejlesztői dokumentáció A szakdolgozat témájában szerepel az ismertetett numerikus módszerek webes megvalósítása. Erre a célra a php nyelvet választottam, ugyanis egy ilyen feladatra szerintem megfelelően alkalmas. Az index.html főoldal 3 frame-re van osztva: felül a szakdolgozat címe, bal oldalon a menü, jobb oldalon pedig a tartalom jelenik meg - minden a menüből megnyitott oldal a ’tartalom’ frame-ben látható. Az index.html -hez kapcsolódik egy .css stílusfájl is, de csak apróságokat tartalmaz, ezeket nem részletezném - az oldal nem a megjelenésre összpontosít, hanem a működésre.
A weblap kezdőoldalán látható a frameset szerkezet: fent a cím, bal oldalon a menü, jobbra a tartalom.
3.1. Bemenő adatok beviteli terve Az optimalizálási feladatok lényege egy adott függvény minimalizálása, vagy feltételek nélkül, vagy valamilyen feltételekkel. Tehát először is egy függvényt kell feldolgozni. A módszerek bemutatásához a hatványfüggvényeket választottam, mivel legtöbbször a tanórák során is azokkal foglalkozunk, valamint ennek csak a deriválás megvalósítása miatt van szerepe, egyébként az optimalizálási eljárások programozását nem befolyásolja ez a szűkítés. A függvények bevitele összeg alakban történik, minden tag felépítése azonos. Meg kell adni először, hogy hány tagból áll a függvény, majd az x vektor dimenzióját. Ebből a program készít egy olyan űrlapot, mely tartalmazza a szükséges beviteli mezőket. 27
A tagok felépítése a következő: az együttható szerepel elöl, majd utána x vektor minden eleméhez meg kell adni egy kitevőt. Például, ha x két dimenziós, a 3 × x2 kifejezést így kell beírni: 3 × x10 x21 . A tagokat minden esetben összeadásjel határolja, így kivonást a negatív együtthatóval lehet elérni. Említettem a deriválást, ez a fenti formulával úgy működik, hogy tagonként végezhetjük el, és amelyik elem szerint deriválunk, annak a kitevőjével meg kell szorozni a tag együtthatóját, majd a kitevőt csökkenteni eggyel. Ez így egészen addig működik, amíg a kitevő 0 nem lesz. Ebből az eljárásból következik az a megkötés, hogy a kitevőknek nemnegatív egész számoknak kell lenniük. Erre viszont nincs ellenőrzés, ugyanis ez a leírás tartalmazza ezt a megkötést, tehát a program használójától elvárható, hogy nem fog ennek ellenére negatív kitevőt használni. Az együtthatóknak nem kötelező egésznek lenniük, lehetnek tizedestörtek is. A program kezdőoldalán, ami az Útmutató menüpont, a bevitellel kapcsolatos fontos pontok le vannak írva, segítő jelleggel. Kitérve az egyes módszerekre, mindegyiknél meg kell adni egy kiindulási vektort, valamint egy ε hibahatárt. Előbbi evidens, hogy miért kell, utóbbi célja pedig az, hogy ha nincs pontos eredmény, akkor egy adott határon belüli közelítés esetén is megáll a program és eredményt ad. Szükség van még a lépések számának maximalizálására annak érdekében, hogy egy adott lépés után mindenképpen megálljon az algoritmus. A SUMT módszernél azonban kis lépésszám kell, mert az ott arra vonatkozik, hogy hányszor növeljük tízszeresére a µ értékét. A Newton módszerhez nem is kell más adat, a kvázi Newton eljárásnál viszont értelemszerűen meg kell adni egy mátrixot is, amely a Hesse mátrixot vagy annak az inverzét fogja helyettesíteni. Hogy pontosan melyiket, arról a használni kívánt módszer dönt, amit a radio típusú input mezővel lehet beállítani. Van azonban emellett egy másik, radio típusú beviteli mező, ami annak eldöntésére szolgál, hogy a Newton, vagy kvázi Newton módszereknél alkalmazzon-e a program vonalmenti minimumkeresést, vagy enélkül, az "alap" algoritmusokat hajtsa végre. Ez a 2.3.3-as alfejezet tartalmához kapcsolódik. A SUMT programnál szükség van a feltételekre is, melyeket a célfüggvényhez hasonló módon kell bevinni, külön megadva az egyenlőséges és az egyenlőtlenséges feltételek számát. µ értéke az általánosságban használt 1-re van állítva, 10-szeres szorzóval. Van egy mező még a kerekítés megadására is, ezzel kapcsolatban szükségesnek érzek egy kis magyarázatot. A nagyon kicsi tizedestörtek nem feltétlenül tudnak pontosan 0 értéket felvenni a float típus pontatlanságai miatt, így pedig nehéz lenne megvizsgálni például azt, hogy a gradiens vektor az x helyen 0 lesz-e. Ezért a megadott számú tizedesjegyre való kerekítéssel oldom meg ezt a problémát, de a többi helyen is a számolás egységesítése miatt, valamint a kiírások során is történik kerekítés. Ez így egy újabb állítható pontosságot biztosít, az ε mellett. Az adatok a form-ból GET metódussal kerülnek feldolgozásra, majd külön tömbben tárolódnak a célfüggvény együtthatói és a kitevői (utóbbiaknak ezen belül több tömb, 28
x indexe szerint). A kitevők sorrendje a következő: $array[$i][$j] az x i-edik eleméhez tartozó kitevő, az összeg j-edik tagjánál. Mivel a dimenziószám, valamint a célfüggvény tagjainak száma változatlan az eljárások alatt, ezért konstansként kerülnek tárolásra.
3.2. Kimenő adatok megjelenítési terve Mind a három módszer külön php fájlban van megírva, azonban sok olyan függvény van, amit több módszer is használ. Ebből kifolyólag egy külön fájlba tettem a közös függvényeket, melyeket így nem kellett minden eljárásnál bemásolni, így átláthatóbb és helytakarékosabb a program. Ez a fuggvenyek.php tartalmazza a kiírással kapcsolatos metódusokat is. Mivel függvényekkel, vektorokkal és mátrixokkal dolgozunk, ezeknek mindenféleképpen kell valamilyen megjelenést adni. A hatványfüggvények kiírására szolgál a legelső, kiir($egy,$ki) nevű programrész. Paraméterként az együtthatók, és a kitevők tömbjét kéri, visszatérési értéke nincs, lényegében a bevitelhez képest olvashatóbb formában adja vissza a függvényeket. Ahol a 0 szerepel, egyszerűsít is, például ha 0 lesz valamelyik együttható, azt a tagot nem írja ki, vagy ha a kitevő 0, a hozzá tartozó xi is láthatatlan marad. A függvény(ek) bevitele után rögtön kiírja a program ezzel a hívással a lap tetejére a célfüggvényt, valamint a feltételeket is, a SUMT esetén. Emellett a gradiens vektor, vagy a Hesse mátrix elemeit is a kiir teszi láthatóvá.
Példa a kiírásra - felül a beviteli mezőknél megadott függvény, alul ugyanez a feldolgozás után megjelenítve. Működése röviden: Megvizsgálja az összes tagot, ha minden együttható 0, akkor kiír egy nullát, ellenkező esetben ha nem nulla az együttható, kiírja az x vektor hozzá tartozó elemeit, a megfelelő kitevőkkel. Ha a kitevők között van nulla, azt az elemet kihagyja. Amennyiben van további nem nulla együttható, akkor + jelet ír elé, majd annak vizsgálatával folytatja. A vektorok és a mátrixok kiírására a vektorkiir($mit) és a matrixkiir($mit) függvények szolgálnak. Az adott dimenzió szerint vektornál egy ciklus, míg mátrixnál két ciklus végigmegy az elemeken, és egy táblázatos formában kiírja azokat. A megjelenést tekintve leginkább a mátrixok determinánsára hasonlít az eredmény, de mivel a témában a determináns nem jelenik meg, ezért nem érthető félre. Ezt a megoldást css segítségével készítettem úgy, hogy a táblázat jobb, és bal oldalán van csak border beállítva. Vektorok azonban megjelenhetnek egyszerű felsorolással is, kerek zárójelben egymás mellé írva, ezt az adott helyen, ha szükség van rá, egy ciklussal oldom meg.
29
Egy dolog van még, ami a megjelenést befolyásolja, mégpedig a módszerek során kapott adatok kiírása. Ahhoz, hogy ne minden egymás alá kerüljön, de jól olvasható legyen egymás mellett, táblázatos megjelenést használok. Minden külön információ saját cellát kap, így elkerülhető a vektorok / mátrixok miatt összefolyt szöveg, és az is, hogy az összefüggő részeket mindig új sor válassza el.
3.3. Programkód Ennyi bevezetés után nézzük mostmár a program fő algoritmusait, és a módszerek megvalósítását. A fuggvenyek.php további metódusait tekintsük át először. Az egyik alapvető függvény a deriválás: derival($egy, $ki, $szerint). function derival($egy, $ki, $szerint){ $eredmeny[1]=$egy; $eredmeny[2]=$ki; $num = count($eredmeny[1]); for ($i=1; $i<=$num; $i++) { if($eredmeny[2][$szerint][$i] != 0) { $eredmeny[1][$i] = $eredmeny[1][$i] * $eredmeny[2][$szerint][$i]; $eredmeny[2][$szerint][$i] = $eredmeny[2][$szerint][$i] - 1; } else $eredmeny[1][$i] = 0; } return $eredmeny; } Ennek a működése egyszerűnek mondható, a paraméterei az együtthatók tömbje, a kitevők (többdimenziós) tömbje és az az index, amelyik szerint deriválni kell. Egy új, több szintű tömb lesz a visszatérési érték, első eleme az együtthatókkal, második a kitevőkkel. A ciklussal az összes tagon végigmegyünk, ha a keresett kitevő nem nulla, akkor az adott tag együtthatóját megszorozzuk a kitevővel, majd a kitevőt csökkentjük eggyel. Ha pedig a kitevő nulla, az együtthatót is nullával tesszük egyenlővé. A gradiens vektor és a Hesse mátrix kiszámításához már csak az előző függvényt kell használni, és ciklikusan meghívni. A gradiens($egy, $ki) létrehozza a gradiens vektort, a hesse($gradiens) függvény pedig a Hesse mátrixot. Érdemes megjegyezni, hogy a hesse függvény paramétere az előzőleg kiszámolt gradiens vektor. Hasonló felépítésű a behelyettesítés megoldása, itt is van egy alap függvény, a helyettesit($egy, $ki, $mit), melynek utolsó paramétere a behelyettesíteni kívánt vektor. Ennél szintén végigmegyünk minden tagon, és elvégezzük a műveleteket, az együtthatót megszorozzuk az aktuális vektor adott hatványaival. Vektorba helyettesítéskor (helyettvekt($mibe, $mit)) egy ciklussal, mátrixnál (helyettmatr($mibe, $mit)) két ciklussal számoljuk ki az összes értéket.
30
Szerepel még egy segédfüggvény - vektorossz($mihez, $mit) - amivel a konkrét értékkel rendelkező vektorokat adhatjuk össze. Ennek például akkor van haszna, ha a Newton módszernél az x vektor következő közelítését számoljuk ki, amikor az aktuális közelítést és az s vektort összeadjuk. A továbbiakban az egyes módszereket részletesen ismertetem.
3.3.1. Newton módszer A modszer/newton.php fájl tartalmazza az algoritmust. Természetesen a közös függvényeket include paranccsal éri el. A form-ról és a szükséges adatokról már esett szó, ezeket nem ismételném. A célfüggvényt a $egyutthatok és a $kitevok tömbökbe tároljuk el, a feldolgozás során. Közben a kiinduló vektort is létrehozzuk, $xvektor néven. Felírjuk a célfüggvény gradiensét - $gradi, majd a Hesse mátrixát - $hes. Ezeket egy táblázatban kiírjuk, az alap adatok után. Itt kezdődik a ciklus, amely az algoritmus minden további részét tartalmazza. 1. Kiszámoljuk a gradiens vektor aktuális értékét az $xvektor-t behelyettesítve $ertekgrad = helyettvekt($gradi, $xvektor); 2. Ha ez a gradiens érték 0, akkor vége a ciklusnak, és kiírjuk a megoldást. 3. Ha a gradiens nem nulla, akkor meg kell oldani a következő egyenletrendszert: H(xk )sk = −∇f (xk ). Ezt kiírja a program az aktuális értékekkel. Az egyenletrendszer megoldása ezzel a függvényhívással történik: $egyrendszer = rendszer($ertekhess, $minuszertekgrad); (a $egyrendszer maga az s vektor értéke) 4. Ha valamiért nincs megoldás, a program befejeződik. Ha van megoldás, akkor a kapott s vektort alapesetben hozzáadjuk az x-hez, ezzel megkapjuk x következő közelítését: $xvektor = vektorossz($xvektor, $egyrendszer); Ha azonban be van állítva a vonalmenti optimalizálás, akkor a s vektort irányként használva (d), a program egy egyváltozós optimalizálást végez, mégpedig Golden section módszerrel. Az eredményként kapott λ érték és a d irány szorzatát adjuk hozzá x-hez. 5. Megnézzük a kilépési feltételt: ha az egymást követő x-ek normája az adott ε alatt van, akkor nincs folytatás, ellenkező esetben folytatjuk a ciklus elejéről. Három függvény kapcsolódik ide, az egyik a kilepfelt($vektor1, $vektor2), ami az x vektorok különbségének első normáját számolja ki, a kilépési feltétel vizsgálatához. A másik az, amivel az egyenletrendszert megoldjuk. A rendszer($bal, $jobb) egy részleges főelemkiválasztásos eljárás programja, mellyel a lineáris egyenletrendszert meg lehet oldani. Ennek az algoritmusa leginkább a lineáris algebra témakörébe tartozik, de röviden a lényege az, hogy a $bal paraméter a Hesse mátrix aktuális értéke, a $jobb pedig a mínusz gradiens vektor, szintén aktualizálva. Ezeket összeillesztem egy új "táblába" ($megold), majd ezen pivotálok, amíg felső háromszög mátrixot kapok, 31
amennyiben ez lehetséges. Ha sikerült, akkor az "eredmény" kommenttel ellátott részben visszafelé haladva kiszámolom az s vektor elemeit.
Példa a Newton módszer egy lépésére. A harmadik függvény pedig a Golden section eljárás megvalósítása: golden($egyutt, $ki, $x, $d, $eps). A paraméterek az erdeti célfüggvény együtthatói és kitevői, az aktuális x és d vektorok, valamint az ε pontossági érték. A függvény működése a következő: először√ is fel kell venni a kiinduló intervallumot $inta, $intb (a és b pontok), majd az α = 5−1 értéket beállítani, és ennek segítségével 2 a belső osztópontokat kiszámolni - $intc és $intd (c és d pontok). Ezután egy ciklus következik, melynek kilépési feltétele az intervallum hosszára vonatkozik: while (($intb-$inta) > 2*$eps) → nagyobb legyen, mint a hibahatár kétszerese. A cikluson belül először az osztópontokkal kell megszorozni a d irányvektort, majd ezt hozzáadni az x-hez, így behelyettesítve a $helyc és $helyd vektorokat a célfüggvénybe, megkapjuk a c és d pontban felvett függvényértéket. Ezeket vizsgálva kell a további lépéseket megtenni a Golden section módszernek megfelelően, az elágaztatás pontosan ezt írja le: Ha f (ck ) > f (dk ) ⇒ [ak+1 , bk+1 ] = [ck , bk ], ck+1 = dk , dk+1 = ak+1 + αLk+1 Ha f (ck ) 5 f (dk ) ⇒ [ak+1 , bk+1 ] = [ak , dk ], ck+1 = ak+1 + (1 − α)Lk+1 , dk+1 = ck . Így egyre kisebb bizonytalansági intervallumokat kapunk, míg végül annak a hossza kisebb lesz, mint a pontossági érték kétszerese - ekkor a végeredményünk egy olyan λ, mely az utolsó intervallum felezőpontja. Megjegyzés: nincs új függvény létrehozva a λ minimalizálására, a φ(λ) = f (xk + λdk ) függvények egyenlőségéből egyértelműen következik, hogy az f függvénybe behelyettesítve az aktuális x és d vektorokat, pontosan egy egyváltozós függvényt kapunk λ-ra, így a φ(λ) külön megalkotására nincs szükség, az eredeti célfüggvény teljesen jól használható. Természetesen nincs különbség a két lehetőség között, a φ(λ) esetén a λ kivételével minden más be van helyettesítve, így az összevonások után kapunk egy ren32
dezett függvényt λ-ra, míg az általam használt módszernél a λ értéke is rögtön bekerül a célfüggvénybe. A λ kiszámolását ez a függvényhívás végzi: $lambda = golden($egyutthatok, $kitevok, $xvektor, $egyrendszer, $epsilon);
3.3.2. Kvázi Newton módszerek A modszer/kvazi.php fájl hajtja végre az algoritmust. Itt is szükség van az include parancsra a függvényekhez. A Newton módszerhez képest annyi változás történik, ami az elméleti összefoglalásban le is van írva, hogy nem a Hesse mátrixszal számolunk, hanem annak, vagy pedig az inverzének egy közelítésével. A beviteli "űrlap" tehát egy mátrix bekérésével lesz több, ebben az esetben. A feldolgozott mátrixot ugyanúgy $hes néven használom, mint az előző módszernél a Hesse mátrixot, így a legegyszerűbb dolgozni vele. A lépések a BFGS eljárás esetén szinte teljesen megegyeznek a Newton módszerrel, az egyetlen különbség, hogy a Hesse mátrixot nem kell meghatározni a gradiens vektorból. A DFP eljárás során viszont a $hes a Hesse mátrix inverzének a közelítését fogja jelölni, és itt a sk = −Dk ∇f (xk ) összefüggést kell megoldani, ahol a D mátrixunk $hes néven szerepel. Viszont ez nem egyenletrendszer, hanem egy vektor szorzása mátrixszal, ami a gyakorlatban egyetlen műveletet jelent. Ezt végzi el a fuggvenyek.php fájlban egy segédfüggvény, az mxv($matrix, $vektor), melynek neve a mátrix szorozva vektor összetételre utal. sk meghatározása után a vonalmenti minimalizálás itt is használható, pontosan ugyanúgy működik, mint a Newton módszernél. A ciklus végén viszont ki kell számolni a helyettesítő mátrix új közelítését is. Ehhez még a következő függvények szerepelnek: vxvt($vektor1, $vektor2) - vektor szorozva vektor transzponáltjával vtxv($vektor1, $vektor2) - vektor transzponáltja szorozva vektorral matrixossz($mihez, $mit) - mátrixok összeadása matrixperszam($matrix, $szam) - mátrix értékeinek leosztása adott számmal A fő függvény pedig az, ami kiszámítja az új közelítő mátrixot: function ujkozelhesse($hesse, $v1, $v2) { $v3 = mxv($hesse, $v2); $szam1 = vtxv($v1, $v2); $mx1 = vxvt($v1, $v1); $per = matrixperszam($mx1, $szam1); $hesse = matrixossz($hesse, $per); $szam2 = -1 * vtxv($v3, $v2); $per2 = matrixperszam(vxvt($v3, $v3),$szam2); $hesse = matrixossz($hesse, $per2); return $hesse; }
33
Ez a néhány függvényhívás pontosan megfelel az elméleti résznél már bemutatott képleteknek. A függvény ráadásul úgy van paraméterezve, hogy bármelyik eljárást választjuk, BFGS-t vagy DFP-t, mindkettőnél működik, mert az alap mátrix ugyanaz lesz, a két vektort pedig csak meg kell cserélni. A következő programrészlet mutatja ezt be: Először is ki kell számolni az y-t: $y = helyettvekt($gradi, $xvektor); $y = vektorossz($y, $minuszertekgrad); $egyrendszer itt is az s vektor értékét jelöli. Ezután pedig a módszereknek megfelelően: if($eljaras == "bfgs") { echo "B<sub>",$lepes+1," =
"; $hes = ujkozelhesse($hes, $egyrendszer, $y); matrixkiir($hes); } Ahogy az fenti kódrészleten is látszik, a két eljárás egy fájlban van megvalósítva, a Newton módszerrel egyező közös részek mellett egy elágaztatás dönti el a szükséges helyeken, hogy melyik metódus esetén milyen lépés következik. A BFGS módszernél lévő ujkozelhesse(\$hes, \$y, \$egyrendszer) függvényhívással bemutatom, hogy az új mátrix felépítése hogyan követi a már ismert képletet: Bk+1 = Bk +
yk ykT (Bk sk )(Bk sk )T − . (Bk sk )T sk ykT sk
ujkozelhesse($hes, $y, $egyrendszer) { $v3 = mxv($hes, $egyrendszer); ⇒ Bk sk $szam1 = vtxv($y, $egyrendszer); ⇒ ykT sk T $mx1 = vxvt($y, $y); ⇒ yk yk y yT $per = matrixperszam($mx1, $szam1); ⇒ ykT sk k k
$hes = matrixossz($hes, $per);
⇒ Bk = Bk +
T yk yk Ts yk k
$szam2 = -1 * vtxv($v3, $egyrendszer); ⇒ −(Bk sk )T sk T k sk )(Bk sk ) $per2 = matrixperszam(vxvt($v3, $v3),$szam2); ⇒ (B−(B T k sk ) sk $hes = matrixossz($hes, $per2);
⇒ Bk = Bk +
(Bk sk )(Bk sk )T −(Bk sk )T sk
A módosított Bk mátrix lesz az új közelítő mátrix, a Bk+1 , mellyel már lehet folytatni a ciklus elejétől a módszert.
34
Példa a DFP kvázi Newton módszer egy lépésére.
3.3.3. SUMT módszer A modszer/sumt.php fájl tartalmazza az algoritmust. A közös függvényeket include paranccsal éri el. A kezdeti két beviteli mező helyett itt öt szerepel, mert meg kell adni az egyenlőtlenséges és az egyenlőséges feltételi függvények számát, valamint ezeknek a maximális tagszámát. Utóbbival kapcsolatban megoldható lett volna, hogy minden feltételnél külön legyen bekérve a tagok száma, de az bonyolította volna a helyzetet, és legtöbbször a tagok száma a feltételeknél közel áll egymáshoz, így egységesen kezelem azokat. Minden ilyen feltételt úgy kell beírni, mint a célfüggvényt. A módszer legfontosabb része, hogy a célfüggvényt minden egyes lépésben kibővítem a megfelelő feltételek "beolvasztásával", de a célfüggvényre később még szükség lesz, ezért nem írom felül, hanem annak adatai helyett segédváltozókat használok. A büntető paraméter a szokásos módon van beállítva: értéke kezdetben 1, és 10szeresére nő minden lépés során. Az adatoknál kért lépésszám fogja megadni, hogy a µ mennyiszer legyen növelve. A programban ezt a következő kód jelzi: $umaxertek = $uertek * pow(10, LEPSZ); - például LEPSZ konstans ha 5, akkor µ legnagyobb értéke 1 ∗ 105 lesz.
35
A büntető paraméter növelése ciklus használatával történik, melynek magjában a feltétel nélküli optimalizálási segédfeladat megoldása szerepel. Ehhez először létre kell hozni a segédfüggvényt, ami a célfüggvény szerepét kapja. Az eredeti célfüggvény együtthatóit és kitevőit egy új változóba másolom, majd az indexek növelésével a tömbhöz fokozatosan csatolom hozzá az éppen szükséges feltételeket. Minden egyes lépésnél először visszaállítom a segédfüggvényt az eredeti állapotára, hogy ne maradjon véletlenül sem több tagja, mint amennyi aktuálisan létrejön. Ennek elérésére az unset parancsot használom a tömbök minden olyan indexű elemére, mely nagyobb, mint az eredeti célfüggvény legnagyobb indexe, vagyis tagjainak száma. Ezután generálom újra a segédfüggvényt. Az egyenlőtlenséges feltételeket, ha bekerülnek, és az egyenlőségeseket pedig minden esetben négyzetre kell emelni. Erre a fvnegyzet($egyutt, $ki) függvény szolgál a fuggvenyek.php fájlból. Működését tekintve nagyon egyszerű: két, egymásba ágyazott ciklus segítségével a függvény minden tagját megszorzom minden taggal. Így például egy négytagú feltételnél a négyzetre emeléskor 16 tagot kapunk, ami egy elég nagy szám, de a program kezeli. Ha elkészült a segédfüggvény, akkor lezajlik a feltétel nélküli optimalizálás. Ennek nem kerül kiírásra minden lépése, csak a végeredmény. A lépésszám itt azonban kötött, konkrétan 100, ha esetleg eléri a program ezt, akkor egy üzenetet ír ki és véget ér az algoritmus. Ha a segédfüggvény minimumát megtalálta a program, akkor egy táblázatba kiírja a büntető paramétert, az x vektort, valamint a segédfüggvény és az eredeti célfüggvény értékét az x helyen. Ezután a µ 10-szeresére nő, és kezdődik újra a segédfeladat, a feltétel nélküli optimalizálás.
Példa SUMT módszerre, 6 lépéssel.
3.3.4. Megjegyzések A program alapvetően a(z előre legyártott) "mintafeladatok" megoldását szolgálja, azokra működik hiba nélkül. Nincsenek ugyanis külön vizsgálatok, például hogy a Hesse mátrix közelítése szimmetrikus-e, azt pedig szintén nem tudni előre, hogy van-e egyáltalán a függvénynek minimumpontja. 36
Egyéb, szabadon választott célfüggvényekre is lehet persze használni, de fontos a végén külön értékelni a kapott eredményt. A kisebb probléma az, hogy többször előfordulhat, hogy utólagos kerekítésre van szükség, például: x1 = −0.99999937 és x2 = 1.00000107 adódik megoldásnak. Persze ez evidens, hogy −1 és 1 lesz (látva, ahogy változnak az x értékek), de ezt tisztázni kell. Hogy a sejtés helyes-e, ki lehet próbálni, ha a kapott végeredményt írjuk be a programnál kezdővektorként. Például lent az látható, hogy a korábbi, kvázi Newton módszerre beillesztett ábrán szereplő feladatnál a kezdővektor helyére a végeredmény triviális kerekítését írtam be (az előző, tizedestörtek helyett a (−1; 1) vektort.
Ellenőrző lépés - egyből megkapjuk a megoldást. A nagyobb probléma az lehet, hogy a függvénynek nem biztos, hogy van minimumpontja. Például az x3 célfüggvényre a Newton módszer előállítja azt a megoldást, hogy x = 0, pedig tudjuk, hogy a függvénynek egyáltalán nincs minimumpontja. De az is előfordulhat, hogy az x vektor értéke egyre növekszik, tehát távolodik az optimális megoldástól. Éppen ezek miatt kell körültekintően elemezni a kapott végeredményt, ha nem pontos. A vonalmenti optimalizálás használatához annyit fűznék hozzá, hogy a kezdő intervallum a [0, 1]-re lett beállítva, de igazából szabadon választható lenne. Nem biztos azonban, hogy a ϕ(λ) függvény a megadott bizonytalansági intervallumon kvázikonvex, ezt figyelembe kell venni. Megemlítenék még egy dolgot, ami a programban történő navigációval kapcsolatos. Ahogy írtam fentebb, lehet úgy ellenőrzést végezni, hogy átírjuk például a kezdővektort, és minden más úgy marad, ahogy volt. Ezt úgy lehet a legegyszerűbben megoldani, hogy a böngészőben a vissza gombot benyomva a beviteli mezőkhöz jutunk, és minden úgy lesz kitöltve, ahogy legutóbb csináltuk. Így könnyen lehet 1-1 dolgot módosítani, akár a kezdővektort, akár a pontosságot, majd újra lefuttatni a programot, vagy az is megoldható, hogy ugyanazt a feladatot megoldjuk vonalmenti keresés használata nélkül, majd visszalépve és átállítva a beviteli mezőt, annak használatával is. Az alkalmazás bemutatása ezzel véget ért, a teljes forráskód megtalálható a mellékelt adathordozón.
37
4. fejezet Összefoglalás A nemlineáris optimalizálás elméletét bemutattam, az alapoktól felépítve. Szó volt a feltétel nélküli optimalizálásról, majd annak szükséges és elégséges feltételeiről. Leírtam ezek után a feladat megoldásának módját, egy példával kísérve. Aztán az ehhez kapcsolódó numerikus módszerek közül kettőt ismertettem. Ahogy említettem ott is, több fajta eljárás létezik, de én a szakdolgozatom témája alapján a Newton típusúakkal foglalkoztam. A Newton módszert, és annak módosulását, a kvázi Newton módszert is összefoglaltam, utóbbinál részletesen kitérve a közelítő mátrixok számolásának lehetőségeire. A másik típusú módszerekről egy rövid említést tettem, mert a vonalmenti minimumkeresés felhasználható a Newton módszer kiegészítéseként, ezt is leírtam. Kellett ehhez egyváltozós optimalizálás is, melyre szintén sok eljárás létezik, én a Golden section-t mutattam be és építettem be a programba. Szintén hoztam példákat is a Newton módszerekre. A feltétel nélküli feladatok mellett a másik fontos része a témának a feltételes optimalizálás. Ezt is az alap definícióktól kezdve építettem fel, a szükséges és elégséges feltételeken át eljutva a KKT feladatig. Ezt megmagyarázva, értelmezve megoldottam egy példát is. Majd a feltételes optimalizálás numerikus módszerei közül soroltam fel néhányat, és részletesen ismertettem egyet, a SUMT metódust, annak is a büntetőfüggvényes változatát. Ehhez kapcsolódó példával zártam az elméleti részt. Ez egy olyan téma, amelynek az elméleti része nem sok újdonságra ad lehetőséget. Éppen ezért van a témának egy másik része is, a webes megjelenítés, ami már igazán önálló munkát és kreativitást igényel. A program megvalósítása a nyelv kiválasztása után a felépítés eldöntésével kezdődött. Én arra jutottam, hogy itt nem a rengeteg féle függvényen van a hangsúly, hanem a módszereken, ezért a függvények csak egy típusára szűkítettem le a feladatot: a hatványfüggvényekre, mert deriválás szempontjából ez az egyik legjobb típus, valamint a tanórákon is rendszerint ilyenekkel foglalkoztunk. Miután kitaláltam, hogy a bevitelük hogyan történjen - tagonkénti megadás, külön tömbökben az együtthatók és a kitevők - meg tudtam írni egy olyan algoritmust, ami elvégzi a deriválást. Innentől már csak követni kellett az elméleti résznél bemutatott módszerek lépéseit, és leprogramozni azokat. Sok segédfüggvényt alkottam, a mátrixokkal és vektorokkal kapcsolatos műveletek, a behelyettesítések, valamint a megjelenítés miatt is. Egy kis
38
gondolkodás után sikerült a saját egyenletrendszer-megoldó függvényemet is létrehoznom. Aztán egy ciklusba tettem a műveleteket, és elkészült a Newton módszer. Ezt módosítva, Hesse mátrix helyett közelítést használva kaptam a kvázi Newton módszereket. Ezeknél persze ügyelni kellett arra, hogy az új közelítő mátrixokat hogyan számoljam ki, ennek is megtaláltam a módját, amit ismertettem a program leírásánál. A SUMT módszer szintén a Newton alapján készült, azt raktam egy újabb, külső ciklusba, valamint megírtam a segédfüggvényt létrehozó sorokat hozzá. Végül a vonalmenti keresést is beépítettem a Newton és kvázi Newton módszerekbe, melyhez megírtam a Golden section módszer algoritmusát. Azt hiszem, mindent, ami a szakdolgozatom témájához tartozott, bemutattam. A program minden eljárását teljesen önerőből, más helyről származó programrészek használata nélkül készítettem, így természetesen ha valahol hibáztam, az is a saját művem, és a félreértés vagy figyelmetlenség számlájára írható. Emellett viszont szerintem jól működnek a módszerek, rendesen végzik a beléjük programozott műveleteket. Persze a végeredményt sokszor még értelmezni kell, mert nem mindig egyértelmű, de a mintafeladatoknál kiszámolja a pontos eredményeket, és ez mutatja, hogy jól működik. A szakdolgozatom itt ér véget, remélem, megfelel az elvárásoknak.
39
Irodalomjegyzék [1] Rapcsák Tamás: Nemlineáris optimalizálás, 2007. (A 2.1-es, Történeti áttekintés résznél lett felhasználva) [2] Dr. Nagy Tamás: Differenciálás, gradiens vektor, Hesse mátrix, láncszabály, implicit függvény tétel, 2012. TÁMOP-4.2.1.B-10/2/KONV-2010-0001 pályázat keretében készült elektronikus tananyag. (Felhasználva a 2.2.1-es alfejezetben) [3] Dr. Nagy Tamás: Feltétel nélküli optimalizálás, 2012. TÁMOP-4.2.1.B-10/2/KONV-2010-0001 pályázat keretében készült elektronikus tananyag. (Felhasználva a 2.2.2-es és 2.2.3-as alfejezetben) [4] Dr. Nagy Tamás: Mátrix definitségének fogalma és tesztek a definitség eldöntésére, 2012. TÁMOP-4.2.1.B-10/2/KONV-2010-0001 pályázat keretében készült elektronikus tananyag. (Felhasználva a 2.2.2-es és 2.2.3-as alfejezetben) [5] Dr. Nagy Tamás: Feltétel nélküli optimalizálás algoritmusai, 2012. TÁMOP-4.2.1.B-10/2/KONV-2010-0001 pályázat keretében készült elektronikus tananyag. (Felhasználva a 2.3-as alfejezetben) [6] Dr. Nagy Tamás: Feltételes optimalizálás, 2012. TÁMOP-4.2.1.B-10/2/KONV-2010-0001 pályázat keretében készült elektronikus tananyag. (Felhasználva a 2.4-es alfejezetben) [7] Dr. Nagy Tamás: Feltételes optimalizálás algoritmusai, 2012. TÁMOP-4.2.1.B-10/2/KONV-2010-0001 pályázat keretében készült elektronikus tananyag. (Felhasználva a 2.5-ös alfejezetben) A definíciók és tételek szó szerinti idézetek (esetleg egy kis átszerkesztéssel), de mivel ezt egyértelműnek érzem, ezért nem használtam idézőjelet. A 2. fejezet többi része tartalmi azonosságokat tartalmaz.
40
Adathordozó használati útmutató A szakdolgozathoz egy CD tartozik mellékletként. 2 mappát tartalmaz: • Program : itt található az elkészült alkalmazás forráskódja. • Szakdolgozat : ez a mappa a szakdolgozatot tartalmazza, LaTeX nyelven A mappák mellett további két dokumentum található a lemezen. Az egyik a szakdolgozat pdf kiterjesztésű változata. A másik pedig a program elindításához és működéséhez kapcsolódó használati útmutató. A program futásához szükséges egy php szerver, egy erre szolgáló alkalmazás hordozható verziójának telepítője van még a lemezen, ez a xampp-portable. A használati útmutatóban le van írva a telepítés menete és a használata.