%% Hetedik gyakorlat
%% Nemlineáris egyenletek: intervallumfelezés, fixpont iteráció
%
% Cél: Az f(x)=0 egyenlet megoldása, ahol f adott, szép tulajdonságú (pl. folytonosan deriválható) függvény, azaz keressük az f függvény gyökeit/zérushelyeit. Ekvivalensen g(x) = b megoldása az f = g-b függvényre alkamazott gyökkereséssel.
% Megjegyzés: A gyökkeresési feladat (ekvivalensen: egyenlet(rendszer) megoldás) hibaérzékenységét nemlineáris függvények esetén is jellemezhetjük egy kondíciószámmal. Egydimenziós esetben, tehát egy f: \mathbb{R} \rightarrow \mathbb{R} függvényre vonatkozólag, egy adott [a,b] intervallumon a zérushelykeresési feladat kondíciószáma legyen
% \kappa = \frac{1}{\inf_{x \in [a,b]}|f'(x)|}.
% Ha ez az érték nagy, akkor a feladat rosszul kondícionált (azaz az f függvény kis megváltoztatásával a megoldás is nagyon megváltozhat).
%% I. Intervallumfelezés
% Tegyük fel, hogy tudjuk, hogy egy adott [a,b] intervallumon a függvényünknek van egy gyöke; például ezt onnan tudhatjuk, hogy a végpontokban különböző előjelű a függvény, pontosabban:f(a) \cdot f(b) < 0. Tegyük fel, hogy ez fennáll és hogy f(a) > 0, f(b)<0. Ez utóbbi nem megszorítás, mert f és -f gyökei ugyanott vannak.
% Alapgondolat: Vizsgáljuk meg a függvény előjelét az \frac{a+b}{2} pontban! Ha itt pozitív, az azt jelenti, hogy a keresett zérushely az \left[a,\; \frac{a+b}{2}\right) intervallumban van, azaz itt érdemes folytatni a keresést. Ugyanakkor ha itt negatív, akkor a zérushely az \left( \frac{a+b}{2},\; b\right] intervallumban található, és itt kell tovább keresnünk. Ezt követően az előbbi gondolatot folytatjuk az új, kisebb intervallumokon, azaz vagy az \frac{a+b}{4}, vagy a \frac{3 (a+b)}{4} pontban vizsgáljuk a előjelet, és így tovább.
% Az algoritmus addig fut, amíg f értéke a vizsgált pontban nulla nem lesz (ez programok esetén a gépi hibák miatt ritkán következik be), vagy amikor a vizsgált intervallum már nagyon kicsi/rövid, de természetesen megadhatunk maximális lépésszámot is, amely után álljon le mindenképpen a program. Határértékben ez a módszer mindenképpen megtalál egy zérushelyet (ha több is van, akkor a zérushelyek közül valamelyiket).
% 1. feladat: Fibonacci egyik cikkében szerepel az alábbi állítás: az
% f(x) = x^3 + 2 x^2 + 10 x - 20
% egyenlet gyöke x = 1.368808107... Ellenőrizzük le az intervallumfelezéses módszerrel, hogy jól számolt-e!
% Megoldás: Legyen a toleranciánk 10^{-10}, és indítsuk az iterációt az [1,2] intervallumon!

format long
tol=10^(-10);

a=1;
b=2;

k=0;

while abs(b-a)>tol
    x = (a+b)/2;
    f = x^3+ 2*x^2+ 10*x-20;
    if f==0
        break; 
    else
        if f>0
            b=x;
        else
            a=x;
        end
    end
    k=k+1;
end
x
k

% Látható, hogy 34 lépés kellett a mószernek, de megtalálta a zérushelyet.
% Használhatunk egy rekurzív megközelítést is:

intervallum_felezes(@(x)(x.^3 + 2*x.^2+ 10*x-20), 1, 2, 1e-10)

% A függvényt is ábrázolhatjuk:

x = linspace(-1,5,101);
f = x.^3+ 2*x.^2+ 10*x-20;
plot(x,f)

%
%% II. Fixpont-iteráció
% Emlék: A Banach-féle fixponttétel következtében ha X Banach-tér és f: X \to X kontrakció 0\leq q < 1 együtthatóval, akkor van x^* fixpontja, amire f(x^*) = x^*; továbbá az x_{n+1} = f(x_n) sorozat tetszőleges X-beli pontból indítva x^*-hoz fog tartani; továbbá
% e_n = x_n - x^* választással látható, hogy
% \|e_{n+1}\| = \| x_{n+1} - x^* \| =  \| f(x_{n}) - f(x^*) \| \leq q \|x_n - x^*\| = q \| e_n\|
% és
% \|e_n\| = \|  (f\circ \ldots \circ f)(x_0) - (f\circ \ldots \circ f)(x^*) \| \leq q^n\|x_0 - x^*\| = q^n\|e_0\|.
% Emlék: Legyen X végesdimenziós valós vektortér, U \subseteq X nyílt halmaz, f: U \to X folytonosan differenciálható és legyen x, h \in Uolyan, hogy az x+h\cdot[0,1] \subseteq U, azaz az x-ből induló x+h végű szakasz része U-nak, akkor
% f(x+h) - f(x) = \int_{0}^1 f'(x+th) \cdot h\, dt.
% Ebből kapjuk a középérték-tétel magasabbdimenziós változatát:
% \|f(x+h) - f(x)\| \leq \left(\sup_{0\leq t \leq 1}\|f'(x+th)\| \right) \|h\|.
% Ha K \subseteq U\subseteq X és Kkompakt, akkor
% M = \sup_{K}\|f'\|
% választással x, y \in K esetén ha az x és y pontok által határolt szakasz K-ban van (pl. ha K konvex), akkor
% \| f(x) - f(y) \| \leq M \|x -y\|.
% Ha f(K) \subseteq K (hívjuk ezt ráképezésnek) és M <1, akkor f|_K:K \to Kkontrakció. Mivel teljes tér zárt részhalmaza teljes, ezért Kteljes, tehát f|_K-ra is alkalmazható a Banach-féle fixponttétel. Érdemes lehet megjegyezni, hogy K általában nem vektortér, azonban az eredeti tér normája által generált metrikával a halmaz teljes metrikus tér, tehát a Banach-féle fixponttétel valóban alkalmazható.
% 2. feladat:
% Oldjuk meg a \cos(x) = x egyenletet a valós számok halmazán. Keressünk alkalmas halmazon alkalmas kontrakciót, majd írjunk kódot amivel addig iteráljunk, míg két szomszédos lépés távolsága nem lesz 10^{-4-nél kisebb.
% Megoldás:
% A baloldal értékkészlete [-1,1], a jobboldal ilyen elemeket csak a [-1,1]-en vesz fel, tehát az egyenletet itt érdemes csak tekinteni. Az oldalak előjelét egyeztetve látjuk, hogy a [0,1] halmazt érdemes néznünk.
% Legyen K=[0,1]. Ekkor f(x) = \cos(x) választással
% 0 \leq f \leq 1 a K halmazon és
% f'(x) = -\sin(x), azaz |f'| < 1 a Kzárt, konvex halmazon - sőt egy kicsit bővebb nyílt halmazon is,
% tehát f: K \to K és kontrakció, szóval alkalmazhatjuk a Banach-féle fixponttételt. Legyen x_0 = 0, hogy ne kelljen válogatni.

x = fixpont_it(0, @(x)(cos(x)), 1e-4, 1000)
cos(x)

% Sajnos az élet nem mindig ad ilyen szép egyenleteket nekünk. Például \cos(6x) = x választással a -6\sin(x) deriváltra már nem adható ilyen egyszerűen jó felső becslés. Általános módszer nincs, de nézzünk meg egy feladatosztályt, amivel el tudunk bírni.
% Monoton (baloldalú) egyenletek
% Tekintsük az g(x) = b feladatot, amire
% g(x) = b \iff f(x) = g(x) - b  = 0,
% sőt c \neq 0 konstans esetén a lánc folytatható: 0 = cf(x) \iff x = x - c f(x) = F_c(x).
% Keressünk egy feltételt, ami mellett valamilyen c választással F_c-re alkalmazható a fixpont-iteráció.
% Tudjuk, hogy F'(x) = I - c f'(x), így ha A = cf'(x), akkor
% \| I- A\|_2^2 = \sup_{\|x\|_2 = 1}x^T(I - A^T)(I - A)x = \sup_{\|x\|_2 = 1}x^T (I - (A+A^T)+ A^TA)x = 1 + \sup_{\|x\|_2 = 1}x^T (-(A+A^T)+ A^TA)x.
% Emiatt
% \|F'(x)\|_2^2 < 1 \iff \sup_{\|x\|_2 = 1}x^T (-(A+A^T)+ A^TA)x < 0,
% átfogalmazva, minden egységhosszú x vektorra
% x^T( -(A+A^T) + A^TA)x  < 0.
% Tegyük fel, hogy A szimmetrikus. EkkorA = Q \Lambda Q^T diagonalizációval és y = Q^Tx helyettesítéssel (ami továbbra is lefedi az összes lehetséges egységvektort)
% -2 y^T \Lambda y + y^T \Lambda^2 y  < 0,
% azaz
% y^T \Lambda(\Lambda - 2I) y < 0.
% Ez nem jelent mást, mint hogy a \Lambda(\Lambda - 2I) mátrix sajátértékei negatívak, szóval minden \lambda sajátértékre \lambda(\lambda - 2) < 0.
% Emélkezzünk, hogy A = cf'(x) volt, tehát implicit feltettük, hogyf'(x) szimmetrikus. Legyen f' szimmetrikus egy K halmaz minden pontjában, továbbá ott a \mu sajátértékei abszolútértékben ne legyenek nagyobbak egy adott M konstansnál. Ekkor tehát minden \mu sajátértékre c\mu(c\mu -2) = c(c\mu^2 - 2\mu) < 0 kellene, hogy teljesüljön.
% Ez a kifejezés c-ben másodfokú, konvex polinom, mely pontosan a gyökei között negatív.
% A gyökök: c_1 = 0  és c_2 = \frac{2\mu}{\mu^2} = \frac{2}{\mu}.
% Láthatjuk, hogy ha van két ellentétes előjelű sajátérték, akkor ez a megközelítés nem működhet (mert c az egyik miatt negatív, a másik miatt pedig pozitív kellene, hogy legyen), tehát tegyük fel, hogy minden sajátérték azonos előjelű (azaz f' a K halmaz minden pontjában (poz. vagy neg.) definit).
% Látható, hogy c_2 a nullához legközelebb akkor van, ha \mu a lehető legnagyobb abszolútértékű (legyen ez az abszolútérték M).
% Ha a sajátértékek negatívak, akkor -\frac{2}{\mu}\leq -\frac{2}{M} < c < 0 jó választás. Ha a sajátérékek pozitívak, akkor 0 < c < \frac{2}{M} \leq \frac{2}{\mu} jó választás.
% Ezekben az esetekben a derivált hosszára vonatkozó felső becslés alapján kontrakciót kapunk az adott K halmazon. Ha F_c(K) \subseteq K és Kmegfelelő (pl. konvex), akkor F_c: K\to Kkontrakció és lesz fixpontja K-ban.
% Összefoglalva: Legyen f egy K zárt, konvex halmazon monoton (azaz f' ott legyen minden pontban (poz. vagy neg.) definit).
% Ha f'szimmetrikus K minden pontjában és a sajátértékeinek abszolútértékére van egy M felső korlát, akkor 0 < c < \frac{2}{M} választással
% ha f poz. monoton K-n, akkor az F_c(x) = x - cf(x),
% ha f neg. monoton K-n, akkor az F_c(x) = x + cf(x)
% függvény kontrakció K-n. Ha ez K-t önmagába is képezi (ráképezés), akkor K \to K kontrakció, tehát K-ban lesz fixpontja és ehhez tartani is tudunk egy K-beli elemből indulva.
% 3. feladat:
% Legyen g(x) = Ax valamilyen poz. def A mátrixra. Tegyük fel, hogy a legnagyobb sajátértéke kisebb, mint egy ismert M konstans. Oldjuk meg az g(x) = b egyenletet fixpont-iterációval. Addig iteráljunk, míg két szomszédos iterált távolsága nem lesz kisebb, mint 10^{-4}.
% Például, vehetjük az alábbi paramétereket:

[Q, R] = qr(magic(3));
A = Q'*diag([1,2,3])*Q
b = ones(3, 1);
M = 3.1;

% de persze lehet mást is.
% Megoldás: Legyen f(x) = Ax - b, ennek a gyökét keressük. f affin, tehát folyt. diffható.
% F_c(x) = x  - cf(x) = x - c(Ax -b), f'(x) = A > 0, tehát f pozit'v monoton. Legyen c pozitív és kisebb, mint \frac{2}{M}. Ekkor f: \mathbb{R}^n\to \mathbb{R}^n kontrakció.
% Ha a fenti érvelést akarjuk használni, akkor mondhatjuk, hogy tetszőlegesen nagy sugarú x^* középpontú zárt gömbben a deriváltjának hossza kisebb mint 1, tehát az iterációval egyet lépve közelebb kerülünk az x^* megoldáshoz; de érvelhetünk közvetlenül az I-cA iterációs mátrix normájával is.

c = 2/M;
x0 = zeros(3, 1);
atol = 1e-4;
max_it = 10^3;

[x, num_steps] = fixpont_it(x0, @(x)(x - c*(A*x -b)), atol, max_it)
A*x

% Érdemes megjegyezni, hogy ez a c választás általában természetesen nem optimális, cserébe működik egy kellően bő egyenletcsaládra. Az ehhez hasonló szimmetrikus poz. def. baloldalú feladatokra használhatjuk a legkisebb és legnagyobb sajátérték számtani közepének reciprokának ellentettjét (ahogy ezt néhány gyakorlattal korábban láttuk).

% c = omega_opt = 2/(lambda_min + lambda_max)
c = 2/(1+3)
[x, num_steps] = fixpont_it(x0, @(x)(x - c*(A*x -b)), atol, max_it)
A*x

% 4. feladat
% Keressük meg az f(x) = x^3 + 2 x^2 + 10 x - 20 egyenlet gyökét, ha valahonnan sejtjük, hogy az a K=[1,2] intervallumban van. Alkalmazzunk fixpont-iterációt és addig iteráljunk, míg a szomszédos iteráltak távolsága 10^{-10} alá nem csökken.
% Megoldás: Mutassuk meg, hogy f monoton K-n!
% Azf'(x) = 3x^2 + 4x + 10 függvény a K halmazon ábrázolva:

xs = linspace(1, 2, 101);
dfxs = 3*xs.^2 + 4*xs + 10;
plot(xs, dfxs)

% tehát sejtjük, hogy itt f' monoton nő. Valóban, deriváltja f''(x) = 6x+4 pozitív ezen az intervallumon. Mivel f' a kezdőpontban pozitív, így f'>0 a K halmazon, tehát f monoton (nő) K-n és így |f'| maximuma a végpontban vétetik fel:
% M = \max_K |f'| = f'(2).
% f monoton nő, tehát valamilyen 0 < c < \frac{2}{M} választás így jó lehet, feltéve, hogy az iteráció nem lép ki a K halmazból.
% Ez utóbbiról meggyőződhetünk analitikusan is, de a számolás helyett elég az is, ha ránézünk a függvény grafikonjára és futtatáskor azt tapasztaljuk, hogy a módszer konvergál.

c = 2/(dfxs(end) + 1)
f = @(x)(x.^3+2*x.^2 +10*x - 20);
plot(xs, xs - c*f(xs))
[x, num_steps] = fixpont_it(1, @(x)(x - c*f(x)), 1e-10, 1000)

%
% Megjegyzés:  Az f(x) = 0 egyenletet természetesen nem csak a korábbiakhoz hasonlóan lehet manipulálni.
% Például, ha f:\mathbb{R} \to \mathbb{R}, f(x) = g(x) x - h(x) alakú, ahol g egy függvény, melynek nincs gyöke K=[a,b]-ban, akkor
% f(x) = 0 \iff x  = \frac{h(x)}{g(x)} = F(x) a Khalmazon.
% Ha F kontrakció K-n és F(K)\subseteq K, akkor alkalmazhatjuk rá a fixpont-tételt a fixpont (azaz f zérushelyének) meghatározásához.
% 5. feladat: Alkalmazzuk ezen megjegyzés ötletét az első feladatban szereplő
% f(x) = x^3 + 2 x^2 + 10 x - 20
% kifejezés [1,2] intervallumban lévő gyökének meghatározására.

%

% Megoldás: Rendezzük át az f(x) = 0 egyenletet az alábbi módon:
% x = \frac{20}{x^2+2x+10}.
% Tehát legyen F(x) = \frac{20}{x^2+2x+10}. Erről kell belátni, hogy az [1,2] intervallumon kontrakció és ráképezés.
% Ezen függvény és deriváltja is monoton az [1,2] intervallumon. Az első feltételhez:
% |F'(x)| = \left| \frac{-20(2x+2)}{(x^2+2x+10)^2} \right| \leq 20 \frac{6}{13^2} = \frac{120}{169}<1,
% ahol felhasználtuk, hogy x \in [1,2]. Hasonlóan, a ráképezésnél:
% 1 \leq \frac{20}{18} \leq \frac{20}{x^2+2x+10}\leq \frac{20}{13}\leq 2,
% tehát valóban ráképezés.

[x, num_steps] = fixpont_it(1, @(x)(20/(x^2+2*x+10)), 1e-10, 1000)

% 6. feladat: Határozzuk meg azon pontot, ahol az f(x)=x^5 - x függvény értéke 2 lesz. Oldjuk meg a feladatot
% intervallumfelezéssel, majd
% fixpont-iterációval!

%

% Megoldás (intervallumfelezés): Azaz most az f(x)=x^5 - x - 2 függvény zérushelyét keressük. Ehhez fel kellene írnunk egy fixpont iterációt. Behelyettesítve néhány értéket kapjuk, hogy f(1) = -2 és f(2) = 28, azaz érdemes az [1,2] intervallumon keresni a gyököt.

x = intervallum_felezes(@(x)(x^5 - x -2), 1, 2, 1e-10)

% Ellenőrzés:

x^5 - x

% Megoldás (fixpont): Azaz most az f(x)=x^5 - x - 2 függvény zérushelyét keressük. Ehhez fel kellene írnunk egy fixpont iterációt. Behelyettesítve néhány értéket kapjuk, hogy f(1) = -2 és f(2) = 28, azaz érdemes az [1,2] intervallumon keresni a gyököt.
% Nulladik ötlet: f monoton nő K=[1,2]-n. M = \max_{x \in [1,2]}|5x^4-1| = 5\cdot2^4 -1 = 79. Próbáljuk ki a c= \frac{2}{M+1} választás mellett F_c-t.

c = 2/(79+1);
[x, num_steps] = fixpont_it(1, @(x)(x-c*(x^5 - x - 2)), 1e-10, 1000)
x^5 - x

% Első ötlet: Átrendezve x = x^5 - 2 adódik, azaz legyen F(x) = x^5 - 2. Ekkor azonban F'(x) = 5 x^4 - 1, ami nem lesz kisebb, mint 1 az [1,2] intervallumon, így nem lesz jó.
% Második ötlet: Fejezzük ki az egyenletben szereplő másik x-et! Ekkor x = \sqrt[5]{x + 2} adódik, azaz legyen F(x) = \sqrt[5]{x + 2}. Ekkor
% |F'(x)| = \left| \frac15 \left( x+2 \right)^{-4/5} \right| = \left| \frac15 \frac{1}{(x+2)^{4/5}} \right| \leq 0.09 < 1,
% azaz kontrakció lesz (használtuk, hogy a függvény monoton az [1,2] intervallumon).
% A ráképezés az alábbi becslésből adódik:
% 1 \leq \sqrt[5] 3 \leq \sqrt[5]{x + 2} \leq \sqrt[5] 4 \leq 2.
% Azaz így mivel ez egy ráképezés és kontrakció, így a fixpont-iteráció konvergálni fog, ha az intervallumból indítjuk.

[x, num_steps] = fixpont_it(1, @(x)((x+2)^(1/5)), 1e-10, 1000)

% Ellenőrzés:

x^5 - x

%% Konvergenciarend
% Definíció: Legyen az x_k \rightarrow x^* konvergens iteráció k-adik elemének hibája e_k (azaz e_k = x_k - x^*). Ekkor azt mondjuk, hogy az iteráció konvergenciarendje p, ha
% p = \lim_{k \rightarrow \infty} \frac{\ln ||e_k||}{\ln ||e_{k-1}||} =\lim_{k \rightarrow \infty} \log_{\|e_{k-1}\|} \|e_k\|,
% tehát határértékben egy hibavektor hosszának durván hányadik hatványa a következő hibavektor hossza.
% Emlék: Fixpont-iterációra a szokásos feltételek mellett
% \|e_{n+1}\| = \| x_{n+1} - x^* \| =  \| F(x_{n}) - F(x^*) \| \leq q \|x_n - x^*\| = q \|e_n\|,
% ahol q az F kontrakció együtthatója. Egyéb feltételek nélkül tehát ez a módszer (legalább) elsőrendű.
% Tétel: Fixpont-iteráció esetén ha F(x_k) \rightarrow x^* és azt is tudjuk, hogy
% F'(x_k) = F''(x_k) = \dots = F^{(p-1)}(x_k)=0 \qquad \text{de} \qquad F^{(p)}(x^*)\neq 0,
% akkor a módszer p-ed rendben konvergál.
% Motiváció:
% e_{k+1} = x_{k+1} - x^* = F(x_k) - F(x^*) = \sum\limits_{j=1}^p \frac1{j!}F^{(j)}(x^*)(x_{k} - x^*)^j + o(\|x_{k} - x^*\|^{p}) =  0+ \frac1{p!}F^{(p)}(x^*) e_k^p + o(\|e_k\|^{p}),
% így
% \|e_{k+1}\| \leq C_p \|e_k\|^p  + o(\|e_k\|^p),
% tehát határértékben \|e_{k+1}\| durván (legalább) p-edik hatványa az \| e_k\| értéknek, azaz p-edrendről beszélhetünk.
% 4. feladat: Tekintsük a következő függvényt: f(x) = \frac{x}{x-1}.
% Hol veszi fel a 2 értéket? Határozzuk meg ezt a pontot kézzel!
% Mit tudunk mondani a konvergenciájának rendjéről, ha fixpont-iterációval közelítenénk?
% Alkalmas intervallum választásával numerikusan ellenőrizzük is a konvergenciarendet!
% Megoldás: Most tehát az f(x) = \frac{x}{x-1} - 2 egyenlet zérushelyét keressük. Viszonylag egyszerűen látszik, hogy
% f(x) = \frac{x}{x-1} - 2 = \frac{-x+2}{x-1} =0 ,
% tehát a megoldás az x=2 pont. A fixpont iteráció alakja ebben az esetben
% F(x)=x + \frac{x}{x-1}-2.
% A tétel alapján a konvegenciarendhez számítsuk ki a deriváltakat az x=2 pontban!
% F'(x) = 1 + \frac{(x-1)-x}{(x-1)^2}, \qquad \text{azaz} \qquad F'(2)=0.
% F''(x)= \frac{2}{(x-1)^3}, \qquad \text{azaz} \qquad F''(2)=2\neq 0.
% Mindezek alapján tehát a módszer másodrendben konvergál.
% Ábrázoljuk a függvényt!

xs = linspace(1.1, 5, 101);
F = @(x)(x + x./(x-1) - 2);

plot(xs, F(xs))

% Látható, hogy ráképezés lesz a [2,3] intervallumon (hiszen F(3)<3); és kontrakció is, mivel a derivált monoton és
% |F'(x)| = \left| 1 - \frac{-1}{(x-1)^2} \right| \leq \frac12 <1.
% Az iteráció:

[x, num_steps] = fixpont_it(3, F, 1e-2, 1000)

% Az iteráció kiegészítve a rendekkel:

% hiba az iteracioval kapott kozelites eseten
e0 = norm(x - 2)
% hiba, ha lepunk meg egyet
e1 = norm(F(x) - 2)
% lehetne tobb lepessel is, de igazabol jo ez is,
% hiszen kellően közel vagyunk már a fixponthoz.
rend_becsles = log(e1)/log(e0)

%
% Szorgalmi: Tekintsük az alábbi függvényt, ahol c \in \mathbb{R}^+:
% f(x) = \frac{1}{x} - \frac{x}{c}
% (a) Keressük a függvény x^*>0 zérushelyét fixpont iteráció segítségével! Írjuk fel alkalmas F(x)  függvényt, és bizonyítsuk be, hogy valóban konvergál minden c>0 esetén! Hová fog tartani?

%

% (b) Mit mondhatunk az iteráció rendjéről c függvényében? Az állításunkat numerikusan is szemléltessük!

%

% Szorgalmi:
% Oldjuk meg az
% \begin{cases}
% x(x^2+y^2) &= 2\\
% y(x^2+y^2) &= -2
% \end{cases}
% egyenletrendszert (síkban haladó) fixpont-iterációval.
% Határozzuk meg a megoldást kézzel.
% Keressünk a megoldást tartalmazó egyszerű alakú, alkalmas Khalmazt, amin a megfelelő függvény monoton. Ezt úgy tudjuk ellenőrizni, hogy kiszámoljuk a deriváltját, majd a derivált sajátértékeit, melyeknek azonos előjelűeknek (és nem nullának) kell lenniük.
% írjunk fel egy alkalmas K esetén egy alkalmas c választás mellett egy F_c függvényt, ami K-t önmagába képezi (ráképezés) és azon kontrakció. A ráképezést nem muszáj ellenőrizni.
% Futtassuk a fixpont-iteráció kódját az adott függvényre egy alkalmas kezdőpontból indítva (ez ne az analitikus megoldás legyen)

%

%
%% Függvénytemető:

intervallum_felezes(@(x)(x.^3 + 2*x.^2+ 10*x-20), 1, 2, 1e-10)
intervallum_felezes(@(x)(x-2), 2-2^15, 2+2^15, 1e-10)
intervallum_felezes(@(x)(x^10), -2^15, 2^15, 1e-10)
[x, num_steps] = fixpont_it(10, @(x)(-x/2), 1e-10, 1000)

%

function r = intervallum_felezes(f, a, b, atol)
    if f(a)*f(b) > 0
        r = 0; return
    end
    
    midpoint = (a+b)/2;
    if abs(a-b) < atol
        r = midpoint; return
    end
    
    % egy nemkisebb abszoluterteku eredmennyel terjunk vissza
    r = max(intervallum_felezes(f, a,        midpoint, atol), ...
            intervallum_felezes(f, midpoint, b,        atol), ...
            'ComparisonMethod', 'abs');
end

%

function [x, num_steps] = fixpont_it(x0, f, atol, max_it)
    x = x0;
    for num_steps = 1:max_it
        dx = -(x - f(x));
        x = x + dx;
        if abs(dx) < atol
            break
        end
    end
end