Entwickler-Ecke

Algorithmen, Optimierung und Assembler - Problem mit RegulaFalsi


ImbaPanda - Di 05.02.08 16:21
Titel: Problem mit RegulaFalsi
Also, ich verwende das oben genannte Verfahren, um Nullstellen einer Gleichung zu bestimmen. Klappt auch soweit ganz gut nur ich habe ein großes Problem. In ca. 80% der Fälle hängt sich das Programm mehr oder weniger auf, weil es einfach zu viele Berechnungen sind. Vielleicht hat jemand eine Idee, wie ich das Ganze einschränken könnte? EIne gewisse Toleranz bei dem Ergebnis wäre ja auch akzeptabel.
Hier mal der Code

Delphi-Quelltext
1:
2:
3:
4:
5:
6:
7:
8:
9:
10:
11:
12:
13:
14:
15:
16:
17:
18:
19:
20:
21:
22:
23:
24:
25:
26:
 
Procedure TForm1.RegulaFalsi(Funktion:String;x1,x2:Real);
var
 Temp,x0: real;
 Expr:Expression;
begin
  Expr:=ParseExpr(Funktion,['x']);
  repeat
    // X0 ist der nächste nähere punkt an der 0-Stelle
    X0 := X1 - EvalExpr(Expr,[X1]) * ( (X2 - X1) / (EvalExpr(Expr,[X2]) - (EvalExpr(Expr,[X1]))) );
    // Wie nah ist X0 an 0 ?
    Temp := EvalExpr(Expr,[X0]);
    if Temp = 0 then
      // gleich 0, fertig !
      Break
    else if Temp > 0 then
      // größer 0 - also näher an X2
      X2 := X0
    else
      // kleiner 0 - also näher an X1
      X1 := X0;
    // stop wenn der abstand zw. x1 und x2 kleiner als 1/1000 ist
  until Abs(x2 - x1) < 0.001;
  //Ergebnis eintragen
     Listbox1.Items.Add(FloatToStr(x0));
end;


alzaimar - Mi 06.02.08 09:15

Verwende ein anderes Verfahren oder analysiere (visualsiere) deine Zwischenschritte. Ich vermute, die Schritte hüpfen wie wild in der Gegend herum, weil die Startwerte nicht geeignet sind.

Welchen math.Ausdruck willst du denn nullen?


Horst_H - Mi 06.02.08 11:57

Hallo,

du setzt voraus, f(x1) < 0 und f(x2) > 0 ist (oder umgekehrt), dies kannst Du ja vorher so bestimmt haben , sieht man aber hier nicht.
Ungleiches Vorzeichen bedeutet, das f(x1)*f(x2)< 0 ist , was sich sehr schnell rechnen laesst.


Delphi-Quelltext
1:
2:
3:
4:
5:
6:
7:
8:
9:
10:
11:
12:
13:
14:
15:
16:
17:
18:
19:
20:
21:
22:
23:
24:
25:
26:
27:
28:
29:
30:
31:
32:
33:
34:
35:
36:
Procedure TForm1.RegulaFalsi(Funktion:String;x1,x2:Real); 
 var 
  x0,f_x0,f_x1,f_x2: double; 
  Expr:Expression; 
 begin 
   Expr:=ParseExpr(Funktion,['x']); 
   //Nicht soviel doppelt rechnen, sondern zwischenspeichern
   //Wer weiss, wie aufwaendig EvalExpr sein kann
   f_x1:= EvalExpr(Expr,[x1]);
   f_x2:= EvalExpr(Expr,[x2]);
   repeat 
     // X0 ist der nächste nähere punkt an der 0-Stelle 
     X0 := X1 -  f_x1* ( (X2 - X1) /(f_x2-f_x2) ); 
     // Wie nah ist X0 an 0 ? 
     f_x0 := EvalExpr(Expr,[X0]); 
     if f_x0 = 0 then 
       // gleich 0, fertig ! 
       Break 
     else 
       if f_x0*f_x2 > 0 then 
         begin
         // gleiches Vorzeichen 
         X2 := X0;
         f_x2 := f_x0; 
         end;
       else 
         begin
         // unterschiedliches Vorzeichen 
         x1 := x0; 
         f_x1 := f_x0;//EDIT war vorher f_x2 := f_x0; Copy&paste
         end;
     // stop wenn der abstand zw. x1 und x2 kleiner als 1/1000 ist 
   until Abs(x2 - x1) < 0.001
   //Ergebnis eintragen 
      Listbox1.Items.Add(FloatToStr(x0)); 
 end;


Ungetestet, da ich gerade PuppyLinux teste...

Gruss Horst


hazard999 - Mi 06.02.08 12:30

Hallo,


Delphi-Quelltext
1:
2:
3:
     f_x0 := EvalExpr(Expr,[X0]); 

     if f_x0 = 0 then


Kommt da ein Integer raus?!?

Wenn nein, ist sicherlich falsch.

Die Abfrage kann schwer greifen, da Fließkomma eigentlich nie 0 sein kann (immer nur beinahe null).


Delphi-Quelltext
1:
Abs(f_x0) < Delta                    


Wäre besser wahrscheinlich besser.

r u

René


ub60 - Mi 06.02.08 12:55

@ImbaPanda:

Soweit ich die Antworten sehe, ist damit Dein Problem gelöst. Wenn nicht und als Tipp für weitere Postings:

ub60


ImbaPanda - Mi 06.02.08 14:38

Also ich hab mich mal drangesetzt und mir was anderes überlegt. Ist relativ einfach und siehe da es funktioniert. Wer Interesse hat hier mal der Code:

Delphi-Quelltext
1:
2:
3:
4:
5:
6:
7:
8:
9:
10:
11:
12:
13:
14:
15:
16:
17:
procedure TForm1.Button1Click(Sender: TObject);
Var
 y1,y2,x,dx,n1:Real;
 Expr:Expression;
begin
x:=1;
dx:=0.00001;
Expr:=ParseExpr(Edit1.Text,['x']);
Repeat
 y1:=EvalExpr(Expr,[x]);
 x:=x + dx;
 y2:=EvalExpr(Expr,[x]);
 If ((y1<0And (y2>0)) Or ((y1>0And (y2<0)) Or (y1=0Or (y2=0)
   Then n1 := x;
Until (n1 = x) or (x > 5);
Listbox1.Items.Add(FLoatToStr(n1));
end;

Es wird simpel überprüft ob zwei direkt aneinander liegende Punkte umgekehrte Vorzeichen besitzen, d.h. ob eine Nullstelle zwischen ihnen ist. Geht recht flott das ganze. Nochmal danke für alle Antworten^^


Horst_H - Mi 06.02.08 16:41

Hallo,

ich habe es mit Lazarus getestest.
Dabei werden einfach nur mal die Startwerte x1,x2 getauscht.
Eine Form mit Button1 und Memo1( schön breit).

Delphi-Quelltext
1:
2:
3:
4:
5:
6:
7:
8:
9:
10:
11:
12:
13:
14:
15:
16:
17:
18:
19:
20:
21:
22:
23:
24:
25:
26:
27:
28:
29:
30:
31:
32:
33:
34:
35:
36:
37:
38:
39:
40:
41:
42:
43:
44:
45:
46:
47:
48:
49:
50:
51:
52:
53:
54:
55:
56:
57:
58:
59:
60:
61:
62:
63:
64:
65:
66:
67:
68:
69:
70:
71:
72:
73:
74:
75:
76:
77:
78:
79:
80:
81:
82:
const
  absDELTA = 1e-6;

type

  { TForm1 }

  TForm1 = class(TForm)
    Button1: TButton;
    Memo1: TMemo;
    procedure Button1Click(Sender: TObject);
  private
    { private declarations }
  public
    { public declarations }
  end

var
  Form1: TForm1; 

implementation

{ TForm1 }

function f(x:double):double;
begin
  result := exp(x)*sin(x);
end;

function RegulaFalsi(x1,x2:double):double;
 var
  x0,f_x0,f_x1,f_x2: double;
//  Expr:Expression;
 begin
     //Expr:=ParseExpr(Funktion,['x']);
   //Nicht soviel doppelt rechnen, sondern zwischenspeichern
   //Wer weiss, wie aufwaendig EvalExpr sein kann
   f_x1:= f(x1);
   f_x2:= f(x2);
   If f_x1*f_x2 > 0 then
     exit;
   
   repeat
     // X0 ist der nächste nähere punkt an der 0-Stelle
     X0 := X1 -  f_x1* ( (X2 - X1) /(f_x2-f_x1) );
     // Wie nah ist X0 an 0 ?
     f_x0 := f(X0);
     FOrm1.Memo1.lines.Add(Format('f(%0.2f) = %0.3f | f(%0.2f) =%0.3f |f(%0.4f) = %e',[x1,f_x1,x2,f_x2,x0,f_x0]));;
     if f_x0 = 0 then
       // gleich 0, fertig !
       Break
     else
       if f_x0*f_x2 > 0 then
         // gleiches Vorzeichen
         begin
         X2 := X0;
         f_x2 := f_x0;
         end
       else
         // unterschiedliches Vorzeichen
         begin
         x1 := x0;
         f_x1 := f_x0;
         end;
     // stop wenn der abstand zw. x1 und x2 kleiner als 1/1000 ist
   until ABS(f_x0) < absDELTA;

   result := x0;
   //Ergebnis eintragen
 end;

procedure TForm1.Button1Click(Sender: TObject);
Var
 x1,x2:double;
begin
  x1 := 2.0;
  x2 := pi+1.0;
  Memo1.lines.Add(Format('%5.2f %5.2f  %e',[x1,x2,RegulaFalsi(x1,x2)]));
  x1 := x2;
  x2 := 2.0;
  Memo1.lines.Add(Format('%5.2f %5.2f  %e',[x1,x2,RegulaFalsi(x1,x2)]));
end;


das scheint zu funktionieren.
Man sieht das Regula falsi sich nur von einer Seite annähert und deshalb langsam konvergiert.
[code}
Memo1
f(2,00) = 6,719 | f(4,14) =-52,931 |f(2,2412) = 7,36921735523644E+000
f(2,24) = 7,369 | f(4,14) =-52,931 |f(2,4735) = 7,34962934533402E+000
f(2,47) = 7,350 | f(4,14) =-52,931 |f(2,6769) = 6,51636966422036E+000
f(2,68) = 6,516 | f(4,14) =-52,931 |f(2,8374) = 5,11314642176185E+000
f(2,84) = 5,113 | f(4,14) =-52,931 |f(2,9523) = 3,60339686804760E+000
f(2,95) = 3,603 | f(4,14) =-52,931 |f(3,0281) = 2,33950391238947E+000
f(3,03) = 2,340 | f(4,14) =-52,931 |f(3,0752) = 1,43598425173683E+000
f(3,08) = 1,436 | f(4,14) =-52,931 |f(3,1034) = 8,50536385191576E-001
f(3,10) = 0,851 | f(4,14) =-52,931 |f(3,1198) = 4,93031846255074E-001
f(3,12) = 0,493 | f(4,14) =-52,931 |f(3,1292) = 2,82204720371674E-001
f(3,13) = 0,282 | f(4,14) =-52,931 |f(3,1346) = 1,60356794200141E-001
f(3,13) = 0,160 | f(4,14) =-52,931 |f(3,1377) = 9,07411042540247E-002
f(3,14) = 0,091 | f(4,14) =-52,931 |f(3,1394) = 5,12266817552414E-002
f(3,14) = 0,051 | f(4,14) =-52,931 |f(3,1403) = 2,88808024246137E-002
f(3,14) = 0,029 | f(4,14) =-52,931 |f(3,1409) = 1,62702992213408E-002
f(3,14) = 0,016 | f(4,14) =-52,931 |f(3,1412) = 9,16215580365822E-003
f(3,14) = 0,009 | f(4,14) =-52,931 |f(3,1414) = 5,15817511144527E-003
f(3,14) = 0,005 | f(4,14) =-52,931 |f(3,1415) = 2,90359531954751E-003
f(3,14) = 0,003 | f(4,14) =-52,931 |f(3,1415) = 1,63434303164253E-003
f(3,14) = 0,002 | f(4,14) =-52,931 |f(3,1416) = 9,19881398648164E-004
f(3,14) = 0,001 | f(4,14) =-52,931 |f(3,1416) = 5,17737999613273E-004
f(3,14) = 0,001 | f(4,14) =-52,931 |f(3,1416) = 2,91395192888954E-004
f(3,14) = 0,000 | f(4,14) =-52,931 |f(3,1416) = 1,64002861180245E-004
f(3,14) = 0,000 | f(4,14) =-52,931 |f(3,1416) = 9,23035935701717E-005
f(3,14) = 0,000 | f(4,14) =-52,931 |f(3,1416) = 5,19499039168691E-005
f(3,14) = 0,000 | f(4,14) =-52,931 |f(3,1416) = 2,92381775942261E-005
f(3,14) = 0,000 | f(4,14) =-52,931 |f(3,1416) = 1,64556681050840E-005
f(3,14) = 0,000 | f(4,14) =-52,931 |f(3,1416) = 9,26148341247565E-006
f(3,14) = 0,000 | f(4,14) =-52,931 |f(3,1416) = 5,21249296329956E-006
f(3,14) = 0,000 | f(4,14) =-52,931 |f(3,1416) = 2,93366386434583E-006
f(3,14) = 0,000 | f(4,14) =-52,931 |f(3,1416) = 1,65110688736126E-006
f(3,14) = 0,000 | f(4,14) =-52,931 |f(3,1416) = 9,29265917176359E-007

2,00 4,14 3,14159261343257E+000

f(4,14) = -52,931 | f(2,00) =6,719 |f(2,2412) = 7,36921735523644E+000
f(4,14) = -52,931 | f(2,24) =7,369 |f(2,4735) = 7,34962934533402E+000
f(4,14) = -52,931 | f(2,47) =7,350 |f(2,6769) = 6,51636966422036E+000
f(4,14) = -52,931 | f(2,68) =6,516 |f(2,8374) = 5,11314642176185E+000
f(4,14) = -52,931 | f(2,84) =5,113 |f(2,9523) = 3,60339686804760E+000
f(4,14) = -52,931 | f(2,95) =3,603 |f(3,0281) = 2,33950391238947E+000
f(4,14) = -52,931 | f(3,03) =2,340 |f(3,0752) = 1,43598425173683E+000
f(4,14) = -52,931 | f(3,08) =1,436 |f(3,1034) = 8,50536385191576E-001
f(4,14) = -52,931 | f(3,10) =0,851 |f(3,1198) = 4,93031846255074E-001
f(4,14) = -52,931 | f(3,12) =0,493 |f(3,1292) = 2,82204720371674E-001
f(4,14) = -52,931 | f(3,13) =0,282 |f(3,1346) = 1,60356794200141E-001
f(4,14) = -52,931 | f(3,13) =0,160 |f(3,1377) = 9,07411042540247E-002
f(4,14) = -52,931 | f(3,14) =0,091 |f(3,1394) = 5,12266817552414E-002
f(4,14) = -52,931 | f(3,14) =0,051 |f(3,1403) = 2,88808024246137E-002
f(4,14) = -52,931 | f(3,14) =0,029 |f(3,1409) = 1,62702992213408E-002
f(4,14) = -52,931 | f(3,14) =0,016 |f(3,1412) = 9,16215580365822E-003
f(4,14) = -52,931 | f(3,14) =0,009 |f(3,1414) = 5,15817511144527E-003
f(4,14) = -52,931 | f(3,14) =0,005 |f(3,1415) = 2,90359531954751E-003
f(4,14) = -52,931 | f(3,14) =0,003 |f(3,1415) = 1,63434303164253E-003
f(4,14) = -52,931 | f(3,14) =0,002 |f(3,1416) = 9,19881398648164E-004
f(4,14) = -52,931 | f(3,14) =0,001 |f(3,1416) = 5,17737999613273E-004
f(4,14) = -52,931 | f(3,14) =0,001 |f(3,1416) = 2,91395192888954E-004
f(4,14) = -52,931 | f(3,14) =0,000 |f(3,1416) = 1,64002861180245E-004
f(4,14) = -52,931 | f(3,14) =0,000 |f(3,1416) = 9,23035935701717E-005
f(4,14) = -52,931 | f(3,14) =0,000 |f(3,1416) = 5,19499039168691E-005
f(4,14) = -52,931 | f(3,14) =0,000 |f(3,1416) = 2,92381775942261E-005
f(4,14) = -52,931 | f(3,14) =0,000 |f(3,1416) = 1,64556681050840E-005
f(4,14) = -52,931 | f(3,14) =0,000 |f(3,1416) = 9,26148341247565E-006
f(4,14) = -52,931 | f(3,14) =0,000 |f(3,1416) = 5,21249296329956E-006
f(4,14) = -52,931 | f(3,14) =0,000 |f(3,1416) = 2,93366386434583E-006
f(4,14) = -52,931 | f(3,14) =0,000 |f(3,1416) = 1,65110688736126E-006
f(4,14) = -52,931 | f(3,14) =0,000 |f(3,1416) = 9,29265917176359E-007

4,14 2,00 3,14159261343257E+000
[/code]

Dein Verfahren eignet sich nur zur Bestimmung von den Startwerten von x1 und x2.
Es sagt nichts über f(x0) aus, nur zwischen x1,x2 eine Nullstelle ist.
Bei Funktionen wie sin(1/x) in der Nähe von 0 wird es aber immer schwierig...

Gruß Horst