Wielowymiarowa Metoda Newtona

Przestrzenie wektorowe, bazy, liniowa niezależność, macierze.... Formy kwadratowe, twierdzenia o klasyfikacji...
janusz47
Użytkownik
Użytkownik
Posty: 7910
Rejestracja: 18 mar 2009, o 16:24
Płeć: Mężczyzna
Podziękował: 30 razy
Pomógł: 1670 razy

Wielowymiarowa Metoda Newtona

Post autor: janusz47 »

Wielowymiarowa Metoda Newtona

Dany jest układ równań nieliniowych z \(\displaystyle{ n - }\) niewiadomymi:

\(\displaystyle{ \begin{cases} f_{1}(x_{1}, x_{2},...,x_{n}) = 0 \\ f_{2}(x_{1}, x_{2},...,x_{n}) = 0 \\.....................\ \ (1)\\ f_{n}(x_{1}, x_{2},...,x_{n}) = 0, \end{cases} }\)

który możemy zapisać w postaci wektorowej

\(\displaystyle{ \vec{f}(\vec{x}) = \vec{0} }\)

gdzie:

\(\displaystyle{ \vec{f}(\vec{x}) = \left[ \begin{matrix} f_{1}(\vec{x})\\ f_{2}(\vec{x}) \\ .....\\ f_{n}(\vec{x})\end{matrix} \right], \ \ \vec{x} = \left[ \begin{matrix} x_{1} \\ x_{2}\\ ...\\ x_{n} \end{matrix} \right]^{T}. }\)

O funkcjach \(\displaystyle{ f_{i}(x_{1},x_{2},...,x_{n}), \ \ i=1,2,...,n }\) zakładamy, że mają ciągłe pochodne pierwszeg rzędu w pewnym obszarze

zawierającym pierwiastek układu równań \(\displaystyle{ (1). }\)

Niech \(\displaystyle{ \vec{x}^{(k)} =(x^{(k)}_{1}, x^{(k)}_{2}, ..., x^{(k)}_{n}) }\) będzie \(\displaystyle{ k - }\) tym przybliżeniem pierwiastka \(\displaystyle{

\vec{\alpha} = (\alpha_{1}, \alpha_{2},...,\alpha_{n}). }\)


Dokładną wartość pierwiastka \(\displaystyle{ \vec{\alpha} }\) można wyrazić w postaci równania:

\(\displaystyle{ \vec{\alpha} =\vec{x}^{(k)} + \vec{\epsilon}^{(k)} .}\)

gdzie:

\(\displaystyle{ \vec{\epsilon}^{(k)} = (\varepsilon^{(k)}_{1}, \varepsilon^{(k)}_{2},...,\varepsilon^{(k)}_{n}) }\)

jest błędem przybliżenia pierwiastka \(\displaystyle{ vec{x}^{(k)} }\) w \(\displaystyle{ k - }\) tej iteracji.

Wobec istnienia ciągłej pochodnej funkcji \(\displaystyle{ \vec{f}(\vec{x}) }\) możemy napisać

\(\displaystyle{ \vec{f}(\vec{\alpha}) \approx \vec{f}^{(k)}(\vec{x}^{(k)}) + \vec{f'}^{(k)}(\vec{x}^{(k)})\cdot \vec{\epsilon}^{(k)} \ \ (1)}\)

Zastępując błąd \(\displaystyle{ \vec{\epsilon}^{(k)} }\) przyrostem \(\displaystyle{ \vec{h}^{(k)} = \vec{x}_{k+1} - \vec{x}_{k} }\) i porównując prawą stronę

równania \(\displaystyle{ (1) }\) do zera, otrzymujemy równanie liniowe

\(\displaystyle{ \vec{f}(\vec{x}^{(k)}) +\cdot (\vec{x}^{(k+1)} - \vec{x}^{(k)})= \vec{0} \ \ (2)}\)

Wzór \(\displaystyle{ (2) }\) jest wzorem rekurencyjnym dla Wielowymiarowej Metody Newtona w postaci macierzowej.

Pochodna \(\displaystyle{ \vec{f'}(\vec{x})}\) jest macierzą Jacobiego.

Zakadając jej nieosobliwość, możemy równanie \(\displaystyle{ (2) }\) napisać w postaci:

\(\displaystyle{ \vec{x}^{(k+1)} = \vec{x}^{(k)} + \vec{f^{-1}}(\vec{x}^{(k)})\cdot \vec{f}(\vec{x}^{(k)}).}\)


Przykład

Rozwiązujemy układ równań:

\(\displaystyle{ \begin{cases} f_{1}(x_{1},x_{2}) = x^3_{1} +3x^2_{2} -21 = 0 \\f_{2}(x_{1},x_{2}) = x^2_{1} +2x_{2} +2 = 0 \end{cases} }\)

Jakobian macierzy układu

\(\displaystyle{ \vec{f'}(x_{1},x_{2}) = \left[ \begin{matrix} 3x^2_{1} & 6x_{2} \\ 2x_{1} & 2 \end{matrix} \right] }\)

Przyjmujemy wektor początkowy \(\displaystyle{ \vec{x}^{(0)} = [1, -1] }\)

Obliczamy dla tego wektora wartość wektor funkcji układu i wartość jej Jakobianu.

\(\displaystyle{ \vec{f}(1,-1) = \left[ \begin{matrix} -17 \\ 1 \end{matrix} \right], \ \ \vec{f'}(1,-1)= \left[ \begin{matrix} 3 & -6 \\ 2 & 2 \end{matrix} \right].}\)

Rozwiązujemy układ równań:

\(\displaystyle{ \left[ \begin{matrix} 3 & -6 \\ 2 & 2 \end{matrix} \right]\cdot \left[ \begin{matrix} h^{(0)}_{1} \\ h^{(0)}_{2} \end{matrix} \right] = - \left[\begin{matrix} -17 \\1 \end{matrix} \right]. }\)

Otrzymujemy

\(\displaystyle{ \vec{h}^{(0)} = \left[ \begin{matrix} 1.555556\\ -2.055560 \end{matrix} \right] }\)

Znajdujemy wektor w następnej iteracji:

\(\displaystyle{ \vec{x}^{(1)} = \left [ \begin{matrix} x^{(0)}_{1} \\ x^{(0)}_{2} \end{matrix}\right] = \left[ \begin{matrix} 1 \\ -1 \end{matrix} \right] + \left[ \begin{matrix} 1.555556\\ -2.055560 \end{matrix} \right] = \left[ \begin{matrix} 2.555556 \\ -3.055556 \end{matrix} \right] }\)

Podobnie znajdujemy pozostałe wektory.

Program w OCTAVE

Kod: Zaznacz cały

function newton_sys(F,JF,xx0,tol,maxit)
  x0=xx0;
  iter=1;
  while(iter<=maxit)
  y=feval(JF,x0)\feval(F,x0);
  xn=x0+y;
  err=max(abs(xn-x0));
  if (err<=tol)
    x=xn;
    fprintf('Metoda Newtona zbieżna po %3.0f  iteracjach do \n',iter)
  x
  return;
else
  x0=xn;
end
iter=iter+1;
endwhile
disp('Metoda Newtona nie jest zbieżna')
x=xn;

function f=f1(X)
  x=X(1);y=X(2);
  f=[x.^3+3*y.^2-21;x.^2+2*y+2];

  function df=df1(X)
    x=X(1);y=X(2);
    df=[3*x.^2,6*y;2*x,2];

>> newton_sys('f1','df1',[1,-1]',10^(8),40)
Metoda Newtona zbieżna po   6 iteracjach do
x =

  1.643038
 -2.349787
 
ODPOWIEDZ