Serwis Edukacyjny
w I-LO w Tarnowie
obrazek

Materiały dla uczniów liceum

  Wyjście       Spis treści       Wstecz       Dalej  

Autor artykułu: mgr Jerzy Wałaszek

©2023 mgr Jerzy Wałaszek
I LO w Tarnowie

W naszym serwisie jest nowszy artykuł o obliczaniu pierwiastków funkcji: "Metody numeryczne".

Metoda siecznych

SPIS TREŚCI
Podrozdziały

Algorytm

Mamy daną funkcję f(x), dwa punkty startowe x1 oraz x2 i przedział [a,b] poszukiwań pierwiastka, do którego należą punkty x1, x2. W przedziale poszukiwań pierwiastka funkcja musi spełniać następujące warunki:

Gdy funkcja f(x) spełnia podane warunki, to w przedziale [a,b] istnieje pierwiastek i możemy go wyszukać metodą siecznych.

Metoda siecznych (ang. secant method) jest ulepszeniem metody regula falsi. W metodzie regula falsi żądamy, aby funkcja f(x) zawsze przyjmowała różne znaki na krańcach przedziału poszukiwań pierwiastka. Różne znaki gwarantują nam istnienie pierwiastka w tym przedziale. Otóż w metodzie siecznych taki wymóg nie jest konieczny (za wyjątkiem pierwszego obiegu). Do wyznaczenia kolejnego przybliżenia pierwiastka funkcji wykorzystujemy dwa poprzednio wyznaczone punkty:

Okazuje się, iż takie podejście do problemu znacznie przyspiesza znajdowanie pierwiastka funkcji. Prześledźmy graficznie tok postępowania:

obrazek Za początkowe punkty x1 i x2 przyjmujemy punkty krańcowe przedziału poszukiwań pierwiastka. Ponieważ funkcja zmienia znak w tym przedziale i jest określona oraz ciągła, mamy gwarancję istnienia pierwiastka.

Z punktów wykresu funkcji dla x1 i x2 prowadzimy sieczną - na rysunku zaznaczoną kolorem czerwonym. Punkt przecięcia się tej siecznej z osią OX daje nam następny punkt x3.

obrazek Drugą sieczną prowadzimy z punktów wykresu funkcji dla x2 oraz x3. Otrzymujemy punkt x4.
obrazek Trzecia sieczna prowadzona jest z punktów wykresu funkcji dla x3 i x4. Otrzymujemy kolejny punkt x5. Już widać, iż leży on bardzo blisko rzeczywistego pierwiastka funkcji.
obrazek Ostatnią sieczną prowadzimy z punktów wykresu funkcji dla x4 i x5. Otrzymany punkt x6 jest już dostatecznie blisko pierwiastka funkcji. Zwróć uwagę, iż kolejne punkty xi zbliżają się do siebie.

Jeśli punkty początkowe x1, x2 zostaną źle dobrane, to algorytm może nie dać rozwiązania. Dlatego stosujemy w nim dodatkowy licznik obiegów pętli wyznaczającej kolejne przybliżenia pierwiastka funkcji. Jeśli po wykonaniu zadanej liczby obiegów (u nas 64) nie zostanie znalezione rozwiązanie, to algorytm przerywa obliczenia z odpowiednim komunikatem. Obliczenia również będą przerwane w przypadku, gdy algorytm stwierdzi, iż sieczna jest równoległa do osi OX. Poprawnym warunkiem zakończenia algorytmu jest:

Na początek:  podrozdziału   strony 

Opis algorytmu

Specyfikacja problemu

Dane wejściowe

f(x) funkcja rzeczywista, której pierwiastek liczymy. Musi być ciągła i określona w przedziale poszukiwań pierwiastka.
x1, x2 punkty krańcowe przedziału poszukiwań pierwiastka funkcji f(x). W dalszej fazie obliczeń punkty te przechowują dwie poprzednie wartości xo - x1 ostatnią, a x2 przedostatnią.    x1, x2R

Dane wyjściowe

xo pierwiastek funkcji f(x).   xoR

Zmienne pomocnicze i funkcje

fo, f1 , f2 – wartości funkcji odpowiednio w punktach xo, x1, x2.    fo, f1 , f2 R
i – licznik obiegów pętli. Obiegi są zliczane wstecz od wartości i = 64.  i ∈ C
εo – określa dokładność porównania z zerem. εo = 0.0000000001
εx – określa dokładność wyznaczania pierwiastka xo. εx = 0.0000000001

Lista kroków

  K01: Czytaj x1 i x2
  K02: f1f(x1);    f2f(x2);
i ← 64
  K03: Dopóki (i > 0) ∧ (| x1 - x2 | > εx):
    wykonuj kroki K04...K08
K04:     Jeśli | f1 - f2 | < εo,
    to pisz "Złe punkty startowe"
    i zakończ
K05:
K06:     Jeśli | fo | < εo
    to idź do kroku K10
K07:     x2x1;    f2f1;
    x1xo;    f1 fo
K08:     i ← i - 1
  K09: Jeśli i = 0,
to pisz "Przekroczony limit obiegów"
i zakończ
  K10: Pisz xo
i zakończ

Schemat blokowy

obrazek

Wykonanie algorytmu siecznych rozpoczynamy od odczytu dwóch punktów x1 oraz x2. Punkty te posłużą do wyznaczenia pierwszej siecznej. Nie muszą one obejmować przedziału, w którym funkcja zmienia znak, jednakże jest to wskazane.

Dla obu punktów obliczamy wartości funkcji umieszczając je w zmiennych odpowiednio f1 i f2. Ponieważ algorytm siecznych może nie dać rozwiązania, stosujemy proste zabezpieczenie w postaci zliczania wykonanych obiegów pętli. Zmienna i pełni rolę licznika tych obiegów. Obiegi zliczane są wstecz. Na początku w zmiennej i określamy maksymalną ilość obiegów, po których algorytm zaprzestanie wykonywać obliczenia pierwiastka.

Rozpoczynamy pętlę obliczającą kolejne przybliżenia pierwiastka funkcji. Pętla kończy się w jednym z trzech przypadków:

  1. Licznik i osiąga zero - wypisujemy odpowiedni komunikat i kończymy algorytm - pierwiastek nie został znaleziony, a xo zawiera zwykle jakieś przypadkowe wartości. Jest to sygnał, iż początkowe punkty x1 i x2 były zbyt odległe od pierwiastka funkcji.
  2. Ostatnie dwa pierwiastki przybliżone x1 i x2 różnią się od siebie o εx lub mniej - w takim przypadku xo jest przybliżonym pierwiastkiem. Wypisujemy xo i kończymy algorytm.
  3. Moduł różnicy wartości funkcji f1 i f2 dla dwóch ostatnich pierwiastków x1 i x2 jest równy zero. Ponieważ różnica ta występuje w mianowniku ułamka we wzorze obliczającym nowy pierwiastek xo, musimy przerwać obliczenia. Wypisujemy odpowiedni komunikat i kończymy algorytm.

Na początku pętli wyznaczamy nowe przybliżenie pierwiastka xo oraz wartość funkcji w tym punkcie fo. Do wyznaczenia xo stosujemy wzór siecznych, który wykorzystuje dwa ostatnio wyliczone przybliżenia pierwiastka x1 i x2.

Jeśli fo jest dostatecznie bliskie zeru, to xo jest przybliżonym pierwiastkiem, wypisujemy xo i kończymy algorytm.

W przeciwnym razie przesuwamy x1 i xo odpowiednio do x2 i x1 wraz z wartościami funkcji w tych punktach i kontynuujemy pętlę. Przesunięcie wartości funkcji f1 i fo zwalnia nas od konieczności ich ponownego wyliczania.

Na początek:  podrozdziału   strony 

Programy

Programy wyznaczają miejsce zerowe funkcji:

Pierwiastków należy poszukiwać w przedziałach [-1,0] i [1,2].

C++
// Program znajduje miejsce zerowe funkcji f(x)
// za pomocą algorytmu siecznych
//---------------------------------------------
// (C)2006 mgr J.Wałaszek       I LO w Tarnowie

#include <iostream>
#include <iomanip>
#include <cmath>
#include <cstdlib>

using namespace std;

const double EPS0 = 0.0000000001; // dokładność porównania z zerem
const double EPSX = 0.0000000001; // dokładność wyznaczenia pierwiastka

// Funkcja, której miejsce zerowe obliczamy
// f(x) = x^3*(x+sin(x^2-1)-1)-1
// <-1,0> i <1,2>
//-----------------------------------------

double f(double x)
{
  return x * x * x * (x + sin(x * x - 1) - 1) - 1;
}

//-----------------------------------------------------
// Program główny
//-----------------------------------------------------

int main()
{

  double x0,x1,x2,f0,f1,f2;
  int    i;

  cout << setprecision(8)     // 8 cyfr po przecinku
       << fixed;              // format stałoprzecinkowy

  cout << "Obliczanie pierwiastka funkcji - metoda siecznych\n"
          "f(x) = x^3*(x+sin(x^2-1)-1)-1\n"
          "-------------------------------------------------\n"
          "(C)2006 mgr Jerzy Walaszek        I LO w Tarnowie\n\n"
          "Podaj dwa wstepne punkty x1 i x2:\n\n";
  cout << "x1 = "; cin >> x1;
  cout << "x2 = "; cin >> x2;
  cout << "\n-------------------------------------------------\n\n"
          "WYNIK:\n\n";
  f1 = f(x1); f2 = f(x2); i = 64;
  while(i && (fabs(x1 - x2) > EPSX))
  {
    if(fabs(f1 - f2) < EPS0)
    {
      cout << "Zle punkty startowe\n";
      i = 0;
      break;
    }
    x0 = x1 - f1 * (x1 - x2) / (f1 - f2);
    f0 = f(x0);
    if(fabs(f0) < EPS0) break;
    x2 = x1; f2 = f1;
    x1 = x0; f1 = f0;
    if(!(--i)) cout << "Przekroczony limit obiegow\n";
  }
  if(i) cout << "x0 = " << setw(15) << x0 << endl;
  cout << "\n---------------------------------------------\n\n";
  system("pause");
  return 0;
}
Pascal
// Program znajduje miejsce zerowe funkcji f(x)
// za pomocą algorytmu siecznych
//---------------------------------------------
// (C)2006 mgr J.Wałaszek       I LO w Tarnowie

program mzf1;

uses math;

const
  EPS0 = 0.0000000001; // dokładność porównania z zerem
  EPSX = 0.0000000001; // dokładność wyznaczenia pierwiastka

// Funkcja, której miejsce zerowe obliczamy
// f(x) = x^3*(x+sin(x^2-1)-1)-1
// <-1,0> i <1,2>
//-----------------------------------------

function f(x : real) : double;
begin
  Result := x * x * x * (x + sin(x * x - 1) - 1) - 1;
end;

//-----------------------------------------------------
// Program główny
//-----------------------------------------------------

var
  x0,x1,x2,f0,f1,f2 : double;
  i                 : integer;

begin
  writeln('Obliczanie pierwiastka funkcji - metoda siecznych');
  writeln('f(x) = x^3*(x+sin(x^2-1)-1)-1');
  writeln('-------------------------------------------------');
  writeln('(C)2006 mgr Jerzy Walaszek        I LO w Tarnowie');
  writeln;
  writeln('Podaj dwa wstepne punkty x1 i x2:');
  writeln;
  write('x1 = '); readln(x1);
  write('x2 = '); readln(x2);
  writeln;
  writeln('-------------------------------------------------');
  writeln('WYNIK:');
  writeln;
  f1 := f(x1); f2 := f(x2); i := 64;
  while (i > 0) and (abs(x1 - x2) > EPSX) do
  begin
    if abs(f1 - f2) < EPS0 then
    begin
      writeln('Zle punkty startowe');
      i := 0;
      break;
    end;
    x0 := x1 - f1 * (x1 - x2) / (f1 - f2);
    f0 := f(x0);
    if abs(f0) < EPS0 then break;
    x2 := x1; f2 := f1;
    x1 := x0; f1 := f0;
    dec(i);
    if i = 0 then writeln('Przekroczony limit obiegow');
  end;
  if i > 0 then writeln('x0 = ',x0:15:8);
  writeln;
  writeln('-------------------------------------------------');
  writeln('Koniec. Nacisnij klawisz Enter...');
  readln;
end.
Basic
' Program znajduje miejsce zerowe funkcji f(x)
' za pomocą algorytmu siecznych
'---------------------------------------------
' (C)2006 mgr J.Wałaszek       I LO w Tarnowie

Declare Function f(ByVal x As Double) As Double

Const EPS0 As Double = 0.0000000001 ' dokładność porównania z zerem
Const EPSX As Double = 0.0000000001 ' dokładność wyznaczenia pierwiastka


'-----------------------------------------------------
' Program główny
'-----------------------------------------------------

Dim As double x0, x1, x2, f0, f1, f2
Dim i As Integer

print "Obliczanie pierwiastka funkcji - metoda siecznych"
print "f(x) = x^3*(x+sin(x^2-1)-1)-1"
print "-------------------------------------------------"
print "(C)2006 mgr Jerzy Walaszek        I LO w Tarnowie"
Print
print "Podaj dwa wstepne punkty x1 i x2:"
Print
Input "x1 = ", x1
Input "x2 = ", x2
Print
print "-------------------------------------------------"
Print
print "WYNIK:"
Print

f1 = f(x1) : f2 = f(x2) : i = 64

While (i > 0) And (Abs(x1 - x2) > EPSX)
   If Abs(f1 - f2) < EPS0 Then
      Print "Złe punkty startowe"
      i = 0
      Exit While
   End If
   x0 = x1 - f1 * (x1 - x2) / (f1 - f2)
   f0 = f(x0)
   If Abs(f0) < EPS0 Then Exit While
   x2 = x1 : f2 = f1
   x1 = x0 : f1 = f0
   i -= 1
   If i = 0 Then Print "Przekroczony limit obiegów"
Wend

If i > 0 Then Print Using "x0 = ######.########"; x0

print
print "-------------------------------------------------"
Print
print "Koniec. Nacisnij klawisz Enter..."

Sleep

End

' Funkcja, której miejsce zerowe obliczamy
' f(x) = x^3*(x+sin(x^2-1)-1)-1
' <-1,0> i <1,2>
'-----------------------------------------

Function f(ByVal x As Double) As Double
   Return x * x * x * (x + Sin(x * x - 1) - 1) - 1
End Function
JavaScript
<html>
  <head>
  </head>
  <body>
    <div align="center">
      <form style="BORDER-RIGHT: #ff9933 1px outset; PADDING-RIGHT: 4px;
                   BORDER-TOP: #ff9933 1px outset; PADDING-LEFT: 4px;
                   PADDING-BOTTOM: 1px; BORDER-LEFT: #ff9933 1px outset;
                   PADDING-TOP: 1px; BORDER-BOTTOM: #ff9933 1px outset;
                   BACKGROUND-COLOR: #ffcc66" name="frmbincode">
        <h3 style="TEXT-ALIGN: center">
          Obliczanie pierwiastka funkcji metodą siecznych
        </h3>
        <p style="TEXT-ALIGN: center">
          <i>f(x) = x<sup>3</sup>(x + sin(x<sup>2</sup> - 1) - 1) - 1</i>
        </p>
        <p style="text-align: center">
          (C)2006 mgr Jerzy Wałaszek I LO w Tarnowie
        </p>
        <hr>
        <p style="text-align: center">
          Wpisz do pól edycyjnych dwa punkty dla pierwszej siecznej
        </p>
        <div align="center">
          <table border="0" cellpadding="8" style="border-collapse: collapse">
            <tr>
              <td>
                x<sub>1</sub> = <input type="text" name="wsp_x1" size="20"
                                       value="1" style="text-align: right">
              </td>
              <td>
                x<sub>2</sub> = <input type="text" name="wsp_x2" size="20"
                                       value="2" style="text-align: right">
              </td>
              <td>
                <input type="button" value="Szukaj pierwiastka" name="B1"
                       onclick="main()">
              </td>
            </tr>
          </table>
        </div>
        <div id="out" align=center>...</div>
      </form>

<script language=javascript>

// Program znajduje miejsce zerowe funkcji f(x)
// za pomocą algorytmu siecznych
//---------------------------------------------
// (C)2006 mgr J.Wałaszek I LO w Tarnowie

var EPS0 = 0.0000000001; // dokładność porównania z zerem
var EPSX = 0.0000000001; // dokładność wyznaczenia pierwiastka

// Funkcja, której miejsce zerowe obliczamy
// f(x) = x^3*(x+sin(x^2-1)-1)-1
// <-1,0> i <1,2>
//-----------------------------------------

function f(x)
{
  return x * x * x * (x + Math.sin(x * x - 1) - 1) - 1;
}

//-----------------------------------------------------
// Program główny
//-----------------------------------------------------

function main()
{

  var x0,x1,x2,f0,f1,f2,i,t;

  x1 = parseFloat(document.frmbincode.wsp_x1.value);
  x2 = parseFloat(document.frmbincode.wsp_x2.value);
  if(isNaN(x1) || isNaN(x2))
    t = "<font color=red><b>Nieprawidłowe dane!</b></font>";
  else 
  {
    f1 = f(x1); f2 = f(x2); i = 64;
    while(i && (Math.abs(x1 - x2) > EPSX))
    {
      if(Math.abs(f1 - f2) < EPS0)
      {
        t = "<font color=red><b>Złe punkty startowe!</b></font>";
        i = 0;
        break;
      }
      x0 = x1 - f1 * (x1 - x2) / (f1 - f2);
      f0 = f(x0);
      if(Math.abs(f0) < EPS0) break;
      x2 = x1; f2 = f1;
      x1 = x0; f1 = f0;
      if(!(--i))
        t = "<font color=red><b>Przekroczony limit obiegów!</b></font>";
    }
    if(i) t = "x0 = " + x0;
  }
  document.getElementById("out").innerHTML = t;
}

</script>
    </div>
  </body>
</html>
Wynik:
Obliczanie pierwiastka funkcji - metoda siecznych
f(x) = x^3*(x+sin(x^2-1)-1)-1
-------------------------------------------------
(C)2006 mgr Jerzy Wałaszek        I LO w Tarnowie

Podaj dwa wstępne punkty x1 i x2:

x1 = 1
x2 = 2

-------------------------------------------------

WYNIK:

x0 = 1,18983299

-------------------------------------------------
Koniec. Naciśnij klawisz Enter...

Tutaj możesz przetestować działanie prezentowanego skryptu. Przykładowa funkcja posiada pierwiastki w przedziałach [-1,0] oraz [1,2].

Obliczanie pierwiastka funkcji metodą siecznych

f(x) = x3(x + sin(x2 – 1) - 1) - 1

(C)2006 mgr Jerzy Wałaszek I LO w Tarnowie


Wpisz do pól edycyjnych dwa punkty dla pierwszej siecznej

x1 = x2 =
...
Na początek:  podrozdziału   strony 

Zespół Przedmiotowy
Chemii-Fizyki-Informatyki

w I Liceum Ogólnokształcącym
im. Kazimierza Brodzińskiego
w Tarnowie
ul. Piłsudskiego 4
©2023 mgr Jerzy Wałaszek

Materiały tylko do użytku dydaktycznego. Ich kopiowanie i powielanie jest dozwolone
pod warunkiem podania źródła oraz niepobierania za to pieniędzy.

Pytania proszę przesyłać na adres email: i-lo@eduinf.waw.pl

Serwis wykorzystuje pliki cookies. Jeśli nie chcesz ich otrzymywać, zablokuj je w swojej przeglądarce.

Informacje dodatkowe.