GNU Octave/Rozwiązanie uwikłanego równania nieliniowego

Rozwiązanie uwikłanego równania nieliniowego

Graficzne rozwiązanie uwikłanego układu równań

Rozwiązać numerycznie uwikłany układ równań:

${\begin{cases}x^{2}+y^{2}=1\\x+y=0\end{cases}}$

Zauważmy, że układ ten opisuje przecięcie okręgu o środku 0 i promieniu 1 z prostą $y=-x$. Istnieją zatem dwa rozwiązania i można je obliczyć dokładnie:

${\begin{cases}x=\pm {\sqrt {2}}/2\approx \pm 0.707107\\y=\mp {\sqrt {2}}/2\approx \mp 0.707107\end{cases}}$

Numerycznie rozwiążemy równanie za pomocą funkcjji fsolve.

Zdefiniujmy pewną funkcję wektorową $h$:

$h{\begin{pmatrix}x\\y\end{pmatrix}}={\begin{pmatrix}x^{2}+y^{2}-1\\x+y\end{pmatrix}}$

Będziemy szukać miejsc zerowych dla h, czyli wektorów ${\begin{pmatrix}x_{0}\\y_{0}\end{pmatrix}}$ takich, że $h{\begin{pmatrix}x_{0}\\y_{0}\end{pmatrix}}={\begin{pmatrix}0\\0\end{pmatrix}}$

Zdefiniujmy funkcję h w pliku h.m:

function z=h(x)
   a=x(1)*x(1)+x(2)*x(2)-1;
   b=x(1)+x(2);
   z=[a,b];
endfunction;

Zdefiniujmy macierz pochodnej h w pliku hprim.m

function [y]=hprim(x)
   y = [2*x(1), 2*x(2); 1, 1 ];
endfunction;

Rozwiązujemy równanie wektorowe funkcją fsolve:

octave:5> fsolve(["h"; "hprim"], [1.0, -1.0])
ans =
  0.70711
 -0.70711
octave:6> fsolve(["h"; "hprim"], [-1.0, 1.0])
ans =
 -0.70711
  0.70711

Widać, że Octave rozwiązał równania poprawnie.