Laboratorium Metod Numerycznych

Ć

w. 5

Temat: Numeryczne rozwiązywanie równań i

układów równań nieliniowych

Ć

wikliński Radosław 163509

Majer Łukasz 163557

1.

Cel ćwiczenia

Praktyczne sprawdzenie wiedzy n/t popularnych metod iteracyjnych rozwiazywania
równań i układów równań nieliniowych. Porównanie przydatności poszczególnych metod
do wyznaczania zer funkcji określonych typów. Prześledzenie związku miedzy rzędem
metody iteracyjnej a szybkością zbieżności ciągu przybliżeń.
2.

Porównanie metody bisekcji, siecznych i Newtona- Raphsona

Funkcja 1: f(x)=arctg(4x)-πx
Pochodna funkcji 1: f’(x)=-π + 4/16x

2

+1

Wartość dokładna: x

*

= 0.25

Punkt startowy dla metody Newtona: x

0

= 0.3

Punkt startowy dla metody siecznych x

0

= 0.3 x

1

= 0.29

Punkt startowy dla metody bisekcji x

a

= 0.225 x

b

= 0.3

M-pliki:
-metoda bisekcji

x=[-2:0.01:6];
y=[(x.^2)-2];
plot(x,y),grid;
b=(0.3);
a=(0.225);
xi=a;

for

i=1:1:20

err(i,:)=[0.25-xi];
hold on
plot(a,0,

'r*'

,b,0,

'b^'

)

if

(((atan(4*b)-pi*b)*(atan(4*xi)-pi*xi)) < 0)

a=xi;
xi=(b+a)/2

else

b=xi;

xi=(b+a)/2

end

end

figure

hold on
plot(xi,0,

'b*'

),grid

figure
zak_x=[1:20];
plot(zak_x,err)

0

0.5

1

1.5

2

2.5

3

3.5

4

-12

-10

-8

-6

-4

-2

0

2

0

2

4

6

8

10

12

14

16

18

20

-0.015

-0.01

-0.005

0

0.005

0.01

0.015

0.02

0.025

-metoda siecznych

x(1)=0.3;
x(2)= 0.29;
c=0.25

for

i=2:1:5

fxi=atan(4*x(i))-pi*x(i);
fxi_1=atan(4*x(i-1))-pi*x(i-1);
x(i+1)=x(i) - (fxi * (x(i) - x(i-1)))/(fxi - fxi_1);
blad(i)=c-x(i+1);

end

plot(x,

'b*'

),grid,xlabel(

'numer przyblizenia'

),ylabel(

'wartosc'

);

plot(blad,

'r-'

),grid,xlabel(

'numer przyblizenia'

),ylabel(

'blad'

);

0

2

4

6

8

10

12

14

0.26

0.28

0.3

0.32

0.34

0.36

numer przyblizenia

1

2

3

4

5

6

7

8

-5

-4.5

-4

-3.5

-3

-2.5

-2

-1.5

-1

-0.5

0

x 10

-3

numer przyblizenia

b

la

d

- metoda Newtona- Raphsona

x(1)=0.3;
c=0.25

for

i=1:1:4

fx=atan(4*x(i))-pi*x(i);
fxx=-pi+(4/(1+16*x(i)^2));
x(i+1)=([x(i)]-[fx/fxx])
blad(i)=c-x(i+1);

end

plot(x,

'b*'

),grid,xlabel(

'numer przyblizenia'

),ylabel(

'wartosc'

)

plot(blad,

'b*'

),grid,xlabel(

'numer przyblizenia'

),ylabel(

'blad'

)

0

1

2

3

4

5

6

7

0.24

0.25

0.26

0.27

0.28

0.29

0.3

numer przyblizenia

w

a

rt

o

s

c

1

2

3

4

5

6

7

8

-6

-5

-4

-3

-2

-1

0

x 10

-3

numer przyblizenia

b

la

d

Funkcja 2: f(x)=2(cos(πx)+ πx)- π
Pochodna funkcji 1: f’(x)= 2 (π-π sin(π x))
Wartość dokładna: x

*

= 0.5

Punkt startowy dla metody Newtona: x

0

= 0.45

Punkt startowy dla metody siecznych x

0

= 0.45 x

1

= 0.46

Punkt startowy dla metody bisekcji x

a

= 0.45 x

b

= 0.6

M-pliki:
-metoda bisekcji

x=[0:0.01:4];
y=[2*(cos(pi*x)+pi*x)-pi];
plot(x,y),grid;
b=(0.45);
a=(0.6);
xi=a;

for

i=1:1:20

err(i,:)=[0.5-xi];
hold on
plot(a,0,

'r*'

,b,0,

'b^'

)

if

(((2*(cos(pi*b)+pi*b)-pi)*(2*(cos(pi*xi)+pi*xi)-pi)) < 0)

a=xi;
xi=(b+a)/2

else

b=xi;

xi=(b+a)/2

end

end

figure

hold on
plot(xi,0,

'b*'

),grid

figure
zak_x=[1:20];
plot(zak_x,err)

0

0.5

1

1.5

2

2.5

3

3.5

4

-5

0

5

10

15

20

25

Wykres błędu:

0

2

4

6

8

10

12

14

16

18

20

-0.1

-0.08

-0.06

-0.04

-0.02

0

0.02

-metoda siecznych

x(1)=0.45;
x(2)= 0.46;
c=0.25

for

i=2:1:8

fxi=2*(cos(pi*x(i))+pi*x(i))-pi;
fxi_1=2*(cos(pi*x(i-1))+pi*x(i-1))-pi;
x(i+1)=x(i) - (fxi * (x(i) - x(i-1)))/(fxi - fxi_1);
blad(i)=c-x(i+1);

end

plot(blad,

'b*'

),grid,xlabel(

'numer przyblizenia'

),ylabel(

'blad'

);

1

2

3

4

5

6

7

8

9

-0.04

-0.02

0

0.02

0.04

0.06

numer przyblizenia

b

la

d

1

2

3

4

5

6

7

8

9

-0.04

-0.02

0

0.02

0.04

0.06

numer przyblizenia

b

la

d

-metoda Newtona-Raphsona

x(1)=0.45;
c=0.5

for

i=1:1:4

fx=2*(cos(pi*x(i))+pi*x(i))-pi;
fxx=2*(pi-pi*sin(pi*x(i)))
x(i+1)=([x(i)]-[fx/fxx])
blad(i)=c-x(i+1);

end

plot(blad,

'b*'

),grid,xlabel(

'numer przyblizenia'

),ylabel(

'blad'

);

0

5

10

15

20

25

0.4

0.45

0.5

0.55

0.6

numer przyblizenia

w

a

rt

o

s

c

1

2

3

4

5

6

7

8

-0.04

-0.03

-0.02

-0.01

0

0.01

0.02

0.03

0.04

0.05

numer przyblizenia

b

la

d

3.

Porównanie pierwiastków wielomianu otrzymanych różnymi metodami:

P=[1 -8 12 50 -250 318 696 -2394 585 9050 -20500 15000];

Pierwiastki otrzymane za pomocą funkcji roots:

5.0000
-3.0000
-2.0000 + 1.0000i
-2.0000 - 1.0000i
1.0000 + 2.0000i
1.0000 - 2.0000i
1.0000 + 2.0000i
1.0000 - 2.0000i
2.0000 + 0.0000i
2.0000 - 0.0000i
2.0000

Pierwiastki otrzymane za pomocą funkcji laguerre.m:

1.9995
1.9995
2.0010
0.9999 - 2.0009i
-2.0000 - 1.0000i
1.0003 + 2.0007i
-2.0000 + 1.0000i
1.0001 - 1.9991i
0.9997 + 1.9993i
-3.0000 - 0.0000i
5.0000 - 0.0000i

Pierwiastki otrzymane za pomocą funkcji bairstow.m:

1.9996 - 0.0004i
1.9996 + 0.0004i
-3.0003
2.0007

Pierwiastki otrzymane za pomocą funkcji lehmer.m:
1.9999 + 0.0000i
4.9999 + 0.0000i
2.1217 + 0.0000i
1.8589 + 0.0000i
0.9655 + 2.0582i
1.0431 + 1.9385i
-1.9972 + 1.0006i
-3.0004 + 0.0013i
-2.0016 - 1.0004i
0.9489 - 2.0761i
1.0611 - 1.9220i

1

2

3

4

5

6

7

8

9

10

11

1.5

2

2.5

3

3.5

4

4.5

5

5.5

Bairstow
Laguerre
Lehmer
Roots

4. Rozwiązanie układu równań za pomocą metody Newtona

Funkcje:

Punkt startowy: (x

0

,y

0

) = (3.3; 0.8)

Wartość poprawna (x

*

,y

*

) = (3 1)

function

funkcje

eps=1*10^-9;

syms

x

y

;

f1=inline(

'x*y-x+y-1'

,

'x'

,

'y'

);

f2=inline(

'x^2+2*y^2-11'

,

'x'

,

'y'

);

F=[f1(x,y); f2(x,y)];

J=jacobian(F,[x y]);

F=inline(F);

J=inline(J);

X=[3.3 0.8 ]';

F(X(1),X(2));

J(X(1),X(2));

X1=X-inv(J(X(1),X(2)))*F(X(1),X(2));

for

k=1:100

X1=X-J(X(1),X(2))\F(X(1),X(2));

if

(norm(X1-X)<eps)

break

end

X=X1;

end

punkt=X1

punkt =

3
1