Matematické Fórum

Nevíte-li si rady s jakýmkoliv matematickým problémem, toto místo je pro vás jako dělané.

Nástěnka
22. 8. 2021 (L) Přecházíme zpět na doménu forum.matweb.cz!
04.11.2016 (Jel.) Čtete, prosím, před vložení dotazu, děkuji!
23.10.2013 (Jel.) Zkuste před zadáním dotazu použít některý z online-nástrojů, konzultovat použití můžete v sekci CAS.

Nejste přihlášen(a). Přihlásit

#1 19. 06. 2015 09:53

rss
Zablokovaný
Příspěvky: 233
Škola: 128.1.1.1
Pozice: root
Reputace:   -42 
 

kubická rovnica, x1 pre výpočet počítačom

Dobrý deň.

Potrebujem vyjadriť reálne x1 kubickej rovnice pre počítačovú simuláciu s presnosťou aspoň 4 platné čísla, viď obrázok. http://s9.postimg.org/6ztvw9bpp/kubikgravitacia.png . Ak by sa dalo tak aj x2 a x3. Nechcem numerickú metódu.

Offline

 

#2 19. 06. 2015 09:53

rss
Zablokovaný
Příspěvky: 233
Škola: 128.1.1.1
Pozice: root
Reputace:   -42 
 

Re: kubická rovnica, x1 pre výpočet počítačom

Neviem či je táto otázka dostatočne ťažká a či som ju zaradil správne. :-)

Offline

 

#3 19. 06. 2015 10:10

Bati
Příspěvky: 2441
Reputace:   191 
 

Re: kubická rovnica, x1 pre výpočet počítačom

Ahoj,
původní Cardanova metoda najde reálný kořen (pokud tedy jde o rovnici s reálnými koeficienty, jinak nemusí reálný kořen existovat). Hned pod tím v tom odkaze je další hezká metoda pro nalezení všech kořenů.

Offline

 

#4 23. 06. 2015 20:04

rss
Zablokovaný
Příspěvky: 233
Škola: 128.1.1.1
Pozice: root
Reputace:   -42 
 

Re: kubická rovnica, x1 pre výpočet počítačom

Bati napsal(a):

Ahoj,
původní Cardanova metoda najde reálný kořen (pokud tedy jde o rovnici s reálnými koeficienty, jinak nemusí reálný kořen existovat). Hned pod tím v tom odkaze je další hezká metoda pro nalezení všech kořenů.

Hneď na úvod som použil Cardanovu metódu a vyjadril X všeobecne, Potrebujem vyjadriť len reálne X1.

Offline

 

#5 23. 06. 2015 23:54

Bati
Příspěvky: 2441
Reputace:   191 
 

Re: kubická rovnica, x1 pre výpočet počítačom

↑ rss:
Jsem přesvědčen, že tvou otázku jsem zodpověděl (↑ Bati:).

Offline

 

#6 24. 06. 2015 11:29

pietro
Příspěvky: 4767
Reputace:   187 
 

Re: kubická rovnica, x1 pre výpočet počítačom

↑ rss:Ahoj, jeden stroj mi vyhodil nasledovné, možno Ti to bude vhod   =>

Offline

 

#7 24. 06. 2015 16:46 — Editoval Rumburak (24. 06. 2015 16:58)

Rumburak
Místo: Praha
Příspěvky: 8691
Reputace:   502 
 

Re: kubická rovnica, x1 pre výpočet počítačom

↑ rss:

Ahoj.

Soudím, že hledáš nějakou numerickou metodu. Kdysi se mi na obdobnou úlohu (také šlo o hledání reálného kořene
kubické rovnice) osvědčila metoda tečen, která konvergovala velmi rychle. Informace o ní jistě na webu najdeš.
Předpokladem jejího využití je mít dostatečnou představu o průběhu funkce $f$ z rovnice $f(x) = 0$  (viz předpoklady
oné metody), což u kubického polynomu není až tak těžké.

Offline

 

#8 24. 06. 2015 17:12

Eratosthenes
Příspěvky: 2801
Reputace:   137 
 

Re: kubická rovnica, x1 pre výpočet počítačom

↑ Rumburak:

jenže ↑ rss: hned na začátku píše "Nechcem numerickú metódu". Škoda. Kdyby ji chtěl, pak takové půlení intervalu je na pět řádků programu a tisícinu vteřiny strojového času. Nechápu proč jsou mu milejší šílenosti typu ↑ pietro:...


Budoucnost patří aluminiu.

Offline

 

#9 24. 06. 2015 22:04

pietro
Příspěvky: 4767
Reputace:   187 
 

Re: kubická rovnica, x1 pre výpočet počítačom

↑ Eratosthenes: Ahoj, smekám pred tvorcami  programov v symbolike .

Offline

 

#10 24. 06. 2015 23:03

rss
Zablokovaný
Příspěvky: 233
Škola: 128.1.1.1
Pozice: root
Reputace:   -42 
 

Re: kubická rovnica, x1 pre výpočet počítačom

↑ pietro: Ďakujem, asi na to budem potrebovať viac času, je to to isté čo vyšlo mne?

Offline

 

#11 25. 06. 2015 08:22

pietro
Příspěvky: 4767
Reputace:   187 
 

Re: kubická rovnica, x1 pre výpočet počítačom

↑ rss:  TenTvoj dáva zaručene vždy reálny koreň ( pre rôzne kombinácie b,d) . Ten "môj reálny" občas vydá ( pri niektorých b,d) aj jeden z komplexných koreňov . Ale keď oba vzorce vydajú reálny koreň, tak sú výsledky u Tvojho a môjho zhodné.

Offline

 

#12 25. 06. 2015 12:34 — Editoval rss (25. 06. 2015 12:36)

rss
Zablokovaný
Příspěvky: 233
Škola: 128.1.1.1
Pozice: root
Reputace:   -42 
 

Re: kubická rovnica, x1 pre výpočet počítačom

pietro napsal(a):

↑ rss:  TenTvoj dáva zaručene vždy reálny koreň ( pre rôzne kombinácie b,d) . Ten "môj reálny" občas vydá ( pri niektorých b,d) aj jeden z komplexných koreňov . Ale keď oba vzorce vydajú reálny koreň, tak sú výsledky u Tvojho a môjho zhodné.

Pri výpočte v OpenOfficeCalc mi vyhadzovalo chybu, záporná odmocnina = komplexné číslo.

Offline

 

#13 25. 06. 2015 18:13

mák
Místo: Vesmír, Galaxie MD
Příspěvky: 898
Reputace:   63 
 

Re: kubická rovnica, x1 pre výpočet počítačom

Zdravím, i když nechceš numerickou metodu, tak ti ji nabídnu. Aspoň se podívej jak je jednoduchá:

Code:

Function CubicFindRealRoot(B#,D#,Z#) As Double
    Dim X As Double 
    X = -9e+99
    While Abs(X-Z)>1e-9
        X = Z
        Z = X-(X^3+B*X^2+D)/(3*X^2+2*B*X)
    Wend
    CubicFindRealRoot = Z
End Function

Předáš koeficienty B a D. Proměnná Z je odhad, čím přesněji jej zadáš, tím se méně cyklů se použije. Vrací nejbližší reálný kořen s přesností na 9 desetinných míst (porovnání While).


LibreOffice Verze: 7.6.6.3, Maxima 5.47.0 (SBCL)

Offline

 

#14 25. 06. 2015 20:02

Eratosthenes
Příspěvky: 2801
Reputace:   137 
 

Re: kubická rovnica, x1 pre výpočet počítačom

↑ mák:

Dva dotazy:

1) Pokud jsem stačil zkontrolovat, funkce řeší rovnici x^3+bx^2+d=0. Co lineární člen? Obecně je přeci  x^3+bx^2+cx+d=0...

2) Jak víš, že to zkonverguje? Zkus dosadit b=2, d=5, z=0.1. Myslím, že se kořene nedobereš...


Budoucnost patří aluminiu.

Offline

 

#15 25. 06. 2015 22:32

mák
Místo: Vesmír, Galaxie MD
Příspěvky: 898
Reputace:   63 
 

Re: kubická rovnica, x1 pre výpočet počítačom

1) Funkce řeší rovnici x^3+B*x^2+D podle zadání uvedeném v prvním příspěvku. Pro obecnou rovnici by byla funkce o něco složitější. Chtěl jsem, aby se výpočet zbytečně nekomplikoval.

2) Pro uvedený příklad b=2, d=5, z=0.1 funkce konverguje k číslu -2,69064744802861, Výsledek (reálný kořen) je podle wxMaximy x=-2.69064744802863. Takže zrovna tady to vyšlo. Je možné, že pro některá čísla to nebude konvergovat. Mám však zkušenost, že Newtonova metoda konverguje většinou velmi dobře a rychle (postačuje pouze několik cyklů) zejména, pokud počáteční nástřel je blízko.


LibreOffice Verze: 7.6.6.3, Maxima 5.47.0 (SBCL)

Offline

 

#16 25. 06. 2015 23:21 — Editoval Eratosthenes (26. 06. 2015 00:00)

Eratosthenes
Příspěvky: 2801
Reputace:   137 
 

Re: kubická rovnica, x1 pre výpočet počítačom

↑ mák:

koukám, že k matematice přistupuješ jako k loterii - když při Člověče nezlob se nechceš jedničku, taky hod většinou vyjde. Newtonova metoda konverguje sice rychle, ale podmínky konvergence jsou dost složité. Ten první příklad jsem napsal špatně - zkus rovnici

x^3+5*x^2+3=0

s "počátečním nástřelem" x=-0.5. A to už vůbec nemluvím o "nástřelu" x=0. Takže do programu bych něco takového určitě nedával.


Budoucnost patří aluminiu.

Offline

 

#17 26. 06. 2015 00:31

mák
Místo: Vesmír, Galaxie MD
Příspěvky: 898
Reputace:   63 
 

Re: kubická rovnica, x1 pre výpočet počítačom

Ten tvůj příklad:

Calc (macro): -5,11467914980952
wxMaxima (reálný kořen): x=-5.114679149809521

Takže to zase vyšlo...

K té loterii:
Pokud musím numericky řešit nějaký výpočet, tak vím v jakém rozmezí hodnot budou vstupní a výstupní hodnoty. Takže vytvořím jednoduchý náhradní výpočet, který mi zaručí, že nástřel se bude lišit od skutečnosti pouze minimálně a Newtonovu metodu použiji pouze pro zpřesnění na požadovaný počet míst. Takto to řeším já a zatím mi to vychází.

Proto jsem řekl, že záleží na počátečním odhadu. To bude muset ale vyřešit někdo jiný.

A nebo zkusit jiný způsob výpočtu.
Já jsem jenom chtěl ukázat, že makro je jednoduché.


LibreOffice Verze: 7.6.6.3, Maxima 5.47.0 (SBCL)

Offline

 

#18 26. 06. 2015 04:26 — Editoval mák (26. 06. 2015 04:27)

mák
Místo: Vesmír, Galaxie MD
Příspěvky: 898
Reputace:   63 
 

Re: kubická rovnica, x1 pre výpočet počítačom

Ještě se k tomu vracím, nedalo mi to a napsal jsem přímou metodu výpočtu všech třech kořenů včetně komplexních.

Code:

Function acos(X#) As Double 
    acos = Atn(-X / Sqr(-X * X + 1)) + 2 * Atn(1) 
End Function

Function QBRT(X As Double) As Double
    QBRT = Abs(X) ^ (1 / 3) * Sgn(X)
End Function

Function CubicRoots(A#, B#, C#, D#) As Variant
    Dim F#,G#,H#,I#,J#,K#,L#,M#,N#,O#,P#,Q#,R#,S#,T#,U#,W#
    Dim Z(2,1) As Variant, Y As Integer 
    For Y=0 To 2
        Z(Y,0) = ""
        Z(Y,1) = ""
    Next Y     
    F = C/A-(B/A)^2/3
    G = (27*A^2*D-9*A*B*C+2*B^3)/(27*A^3)
    H = G^2/4+F^3/27
    If H<0 Then ' Tři různé reálné kořeny
        I = sqr(G^2/4-H)
        J = QBRT(I)
        K = acos(-G/(2*I))
        M = cos(K/3)
        N = sqr(3)*sin(K/3)
        P = -B/(3*A)
        Z(0,0) = 2*J*M+P
        Z(1,0) = P-J*(M+N)
        Z(2,0) = P-J*(M-N)
    ElseIf H>0 Then ' Jeden reálný kořen a dva komplexní
        R = -G/2
        T = sqr(H)
        S = QBRT(R+T)
        U = QBRT(R-T)
        L = S+U
        P = -B/(3*A)
        Z(0,0) = +L  +P
        Z(1,0) = -L/2+P
        Z(2,0) = Z(1,0)
        Z(1,1) = +(S-U)*sqr(3)/2
        Z(2,1) = -Z(1,1)
    Else ' Tři identické reálné kořeny
        W = QBRT(-(D/A))
        Z(0,0) = W
        Z(1,0) = W
        Z(2,0) = W    
    EndIf
    CubicRoots = Z        
End Function

Funkce vrací pole 3x2 (levý sloupeček reálná hodnota, pravý imaginární) tří kořenů. Funkce se zadává maticově, takže klasicky trojhmatem [Control]+[Shift]+[Return].
Počítá obecnou rovnici A*X^3+B*X^2+C*X+D = 0.


LibreOffice Verze: 7.6.6.3, Maxima 5.47.0 (SBCL)

Offline

 

#19 26. 06. 2015 11:20

mák
Místo: Vesmír, Galaxie MD
Příspěvky: 898
Reputace:   63 
 

Re: kubická rovnica, x1 pre výpočet počítačom

Demonstrační program ke stažení zde.


LibreOffice Verze: 7.6.6.3, Maxima 5.47.0 (SBCL)

Offline

 

#20 02. 07. 2015 13:31

Honzc
Příspěvky: 4599
Reputace:   244 
 

Re: kubická rovnica, x1 pre výpočet počítačom

↑ rss:
V Pascalu
Výstup pole stringů (kořenů)

Code:

//řeší rovnici: ax^3+bx^2+cx+d=0
//do uses přidat knihovnu Math

type 
  TAlgr=array of string;

//pres - zadá se na kolik platných číslic chceme výsledek (číslo, např.8)

function Algr3(var a,b,c,d: Extended; pres: Integer): TAlgr;
var p,q,dis,y1,y2,y3,yr,yi,x1,x2,x3,fi,r1: Extended;
  rc,ic1,ic2,sy,sx,pr: string;
begin
  SetLength(Result,3);
  pr := '%.'+IntToStr(pres)+'g';
  if a=0 then
  begin
    Result[0] := 'Rovnice není 3.stupně';
    Result[1] := '';
    Result[2] := '';
    MessageBeep($FFFFFFFF);
    Exit;
  end;
    p := (3*a*c-b*b)/(3*a*a);
  if abs(p)<1E-14 then
    p := 0;
  q := (2*b*b*b)/(27*a*a*a)-(b*c)/(3*a*a)+(d/a);
  if abs(q)<1E-14 then
    q := 0;
  dis := -Power(p/3,3)-Power(q/2,2);
  if abs(dis)<1E-15 then
    dis := 0;
  if q>=0 then
    r1 := Sqrt(abs(p/3))
  else
    r1 := -Sqrt(abs(p/3));
  if p<0 then
  begin
    if dis>=0 then
    begin
      try
        if abs(q/(2*Power(r1,3))-1)<1E-14 then
          fi := 0
        else
          fi := ArcCos(q/(2*Power(r1,3)));
        y1 := -2*r1*Cos(fi/3);
        if abs(y1)<1E-14 then
          y1 := 0;
        y2 := 2*r1*Cos(Pi/3-fi/3);
        if abs(y2)<1E-14 then
          y2 := 0;
        y3 := 2*r1*Cos(Pi/3+fi/3);
        if abs(y3)<1E-14 then
          y3 := 0;
        sx := Format('%.15g', [b/(3*a)]);
        sy := Format('%.15g', [y1]);
        x1 := StrToFloat(sy)-StrToFloat(sx);
        if abs(x1)<1E-14 then
          x1 := 0;
        sy := Format('%.15g', [y2]);
        x2 := StrToFloat(sy)-StrToFloat(sx);
        if abs(x2)<1E-14 then
          x2 := 0;
        sy := Format('%.15g', [y3]);
        x3 := StrToFloat(sy)-StrToFloat(sx);
        if abs(x3)<1E-14 then
          x3 := 0;
        Result[0] := Format(pr, [x1]);
        Result[1] := Format(pr, [x2]);
        Result[2] := Format(pr, [x3]);
      except
        on EInvalidOp do
        begin
          Result[0] := 'Přetečená paměť!';
          Result[1] := '';
          Result[2] := '';
        end;
      end;
    end
    else
    begin
      try
        if abs(q/(2*Power(r1,3))-1)<1E-14 then
          fi := 0
        else
          fi := ArcCosh(q/(2*Power(r1,3)));
        y1 := -2*r1*Cosh(fi/3);
        if abs(y1)<1E-14 then
          y1 := 0;
        yr := r1*Cosh(fi/3);
        if abs(yr)<1E-14 then
          yr := 0;
        yi := Sqrt(3)*r1*Sinh(fi/3);
        if abs(yi)<1E-14 then
          yi := 0;
        sx := Format('%.15g', [b/(3*a)]);
        sy := Format('%.15g', [y1]);
        x1 := StrToFloat(sy)-StrToFloat(sx);
        if abs(x1)<1E-14 then
          x1 := 0;
        sy := Format('%.15g', [yr]);
        x2 := StrToFloat(sy)-StrToFloat(sx);
        if abs(x2)<1E-14 then
          x2 := 0;
        if (x2=0) then
        begin
          if (yi=0) then
          begin
            rc := Format(pr, [x2]);
            ic1 := '';
            ic2 := '';
          end
          else
          begin
            rc := '';
            ic1 := '- i '+Format(pr, [abs(yi)]);
            ic2 := 'i '+Format(pr, [abs(yi)]);
          end;
        end
        else
        begin
          rc := Format(pr, [x2]);
          if (yi=0) then
          begin
            ic1 := '';
            ic2 := '';
          end
          else
          begin
            ic1 := ' - i '+Format(pr, [abs(yi)]);
            ic2 := ' + i '+Format(pr, [abs(yi)]);
          end;
        end;
        Result[0] := Format(pr, [x1]);
        Result[1] := rc+ic1;
        Result[2] := rc+ic2;
      except
        on EInvalidOp do
        begin
          Result[0] := 'Přetečená paměť!';
          Result[1] := '';
          Result[2] := '';
        end;
      end;
    end;
  end
  else
  begin
    if p<>0 then
    begin
      try
        if abs(q/(2*Power(r1,3)))<1E-15 then
          fi := 0
        else
          fi := ArcSinh(q/(2*Power(r1,3)));
        y1 := -2*r1*Sinh(fi/3);
        if abs(y1)<1E-14 then
          y1 := 0;
        yr := r1*Sinh(fi/3);
        if abs(yr)<1E-14 then
          yr := 0;
        yi := Sqrt(3)*r1*Cosh(fi/3);
        if abs(yi)<1E-14 then
          yi := 0;
        sx := Format('%.15g', [b/(3*a)]);
        sy := Format('%.15g', [y1]);
        x1 := StrToFloat(sy)-StrToFloat(sx);
        if abs(x1)<1E-14 then
          x1 := 0;
        sy := Format('%.15g', [yr]);
        x2 := StrToFloat(sy)-StrToFloat(sx);
        if abs(x2)<1E-14 then
          x2 := 0;
        if (x2=0) then
        begin
          if (yi=0) then
          begin
            rc := Format(pr, [x2]);
            ic1 := '';
            ic2 := '';
          end
          else
          begin
            rc := '';
            ic1 := '- i '+Format(pr, [abs(yi)]);
            ic2 := 'i '+Format(pr, [abs(yi)]);
          end;
        end
        else
        begin
          rc := Format(pr, [x2]);
          if (yi=0) then
          begin
            ic1 := '';
            ic2 := '';
          end
          else
          begin
            ic1 := ' - i '+Format(pr, [abs(yi)]);
            ic2 := ' + i '+Format(pr, [abs(yi)]);
          end;
        end;
        Result[0] := Format(pr, [x1]);
        Result[1] := rc+ic1;
        Result[2] := rc+ic2;
      except
        on EInvalidOp do
        begin
          Result[0] := 'Přetečená paměť!';
          Result[1] := '';
          Result[2] := '';
        end;
      end;
    end
    else
    begin
      if ((a>0)and(d>0))or((a<0)and(d<0)) then
        r1 := -1
      else
        r1 := 1;
      y1 := r1*Power(abs(q),1/3);
      if abs(y1)<1E-14 then
        y1 := 0;
      x1 := y1-b/(3*a);
      if (q<>0) then
      begin
        yr := -y1/2-b/(3*a);
        yi := y1*Sqrt(3)/2;
        ic1 := ' - i '+Format(pr, [abs(yi)]);
        ic2 := ' + i '+Format(pr, [abs(yi)]);
        rc := Format(pr, [yr]);
        Result[0] := Format(pr, [x1]);
        Result[1] := rc+ic1;
        Result[2] := rc+ic2;
      end
      else
      begin
        Result[0] := Format(pr, [x1]);
        Result[1] := Format(pr, [x1]);
        Result[2] := Format(pr, [x1]);
      end;
    end;
  end;
end;

Offline

 

Zápatí

Powered by PunBB
© Copyright 2002–2005 Rickard Andersson