Ansatzkoeffizienten berechnen

[jan001]

Für n=20 Punkte P(x,y) möchte ich Ansatzfunktionen finden.
Bei Funktion y = a1 + a2*x + a3*x*x + a4*x*x*x + a5*x*x*x*x habe ich eine Lösung zur Berechnung Koeffizienten a1 bis a5 gefunden (Polynom 4. Ordnung).
Für Funktion y = a1 + a2*x + a3*x*x + a4*x*x*x + a5*sqrt(x) hilft diese Methode nicht. Kennt jemand eine Lösung?

Code:

procedure interp_5(n : integer; p : rp; var a : array of double);
var i : integer;
    c1, c2, c3, c4, c5 : double;
    x0, xi, xi2, xi3, xi4, xi5, xi6, xi7, xi8,
    yi, xiyi, xi2yi, xi3yi, xi4yi : double;
begin
x0:= 0; xi:= 0; xi2:= 0; xi3:= 0; xi4:= 0; xi5:= 0; xi6:= 0; xi7:= 0; xi8:= 0;
yi:= 0; xiyi:= 0; xi2yi:= 0; xi3yi:= 0; xi4yi:= 0;

for i:= 1 to n do begin
  x0:= x0 + 1;
  xi:= xi + p.x[i];
  xi2:= xi2 + p.x[i]*p.x[i];
  xi3:= xi3 + p.x[i]*p.x[i]*p.x[i];
  xi4:= xi4 + p.x[i]*p.x[i]*p.x[i]*p.x[i];
  xi5:= xi5 + p.x[i]*p.x[i]*p.x[i]*p.x[i]*p.x[i];
  xi6:= xi6 + p.x[i]*p.x[i]*p.x[i]*p.x[i]*p.x[i]*p.x[i];
  xi7:= xi7 + p.x[i]*p.x[i]*p.x[i]*p.x[i]*p.x[i]*p.x[i]*p.x[i];
  xi8:= xi8 + p.x[i]*p.x[i]*p.x[i]*p.x[i]*p.x[i]*p.x[i]*p.x[i]*p.x[i];
  yi:= yi + p.y[i];
  xiyi:= xiyi + p.x[i]*p.y[i];
  xi2yi:= xi2yi + p.x[i]*p.x[i]*p.y[i];
  xi3yi:= xi3yi + p.x[i]*p.x[i]*p.x[i]*p.y[i];
  xi4yi:= xi4yi + p.x[i]*p.x[i]*p.x[i]*p.x[i]*p.y[i];
end;

Loesung_LGS_5(x0,xi,xi2,xi3,xi4,xi5,xi6,xi7,xi8, yi,xiyi,xi2yi,xi3yi,xi4yi, c1, c2, c3, c4, c5);
a[1]:= c1; a[2]:= c2; a[3]:= c3; a[4]:= c4; a[5]:= c5;
end;

procedure Loesung_LGS_5(x0,xi,xi2,xi3,xi4,xi5,xi6,xi7,xi8, yi,xiyi,xi2yi,xi3yi,xi4yi : double; var c1, c2, c3, c4, c5 : double);
const n = 5;
var
  i0 : integer;
  x : array[0..n-1] of extended;             //Die Lösungen x0, x1, x2, ...
  koeff : array[0..n-1,0..n] of extended;    //Die Koeffizienten

  procedure VereinfacheRest(i: integer);
  var
    zeile, spalte: integer;
    d: extended;
  begin
    for zeile:= 0 to n - 1 do if zeile <> i then begin
      d:= koeff[zeile,i]/koeff[i,i];
      if d <> 0 then begin
        for spalte:= 0 to n do if spalte <> i then
          koeff[zeile,spalte]:= koeff[zeile,spalte] - d*koeff[i,spalte]
        else koeff[zeile,i]:= 0;
      end;
    end;
  end;

begin
koeff[0,0]:= x0;  koeff[0,1]:= xi;  koeff[0,2]:= xi2; koeff[0,3]:= xi3; koeff[0,4]:= xi4; koeff[0,5]:= yi;
koeff[1,0]:= xi;  koeff[1,1]:= xi2; koeff[1,2]:= xi3; koeff[1,3]:= xi4; koeff[1,4]:= xi5; koeff[1,5]:= xiyi;
koeff[2,0]:= xi2; koeff[2,1]:= xi3; koeff[2,2]:= xi4; koeff[2,3]:= xi5; koeff[2,4]:= xi6; koeff[2,5]:= xi2yi;
koeff[3,0]:= xi3; koeff[3,1]:= xi4; koeff[3,2]:= xi5; koeff[3,3]:= xi6; koeff[3,4]:= xi7; koeff[3,5]:= xi3yi;
koeff[4,0]:= xi4; koeff[4,1]:= xi5; koeff[4,2]:= xi6; koeff[4,3]:= xi7; koeff[4,4]:= xi8; koeff[4,5]:= xi4yi;

for i0:= 0 to n - 1 do VereinfacheRest(i0);
for i0:= 0 to n - 1 do x[i0]:= koeff[i0,n]/koeff[i0,i0];

c1:= x[0]; c2:= x[1]; c3:= x[2]; c4:= x[3]; c5:= x[4];
end;

[Furtbichler]

Oh, es gibt ja diverse numerische Näherungsverfahren. Mir fällt sofort Regula Falsi oder Newton ein. Allerdings ist das bei Polynomen höherer Ordnung so eine Sache, denn Du kannst keine beliebigen Startwerte nehmen. Aber Wenn Du einigermaßen weißt, in welchem Rahmen sich die Lösung bewegt, kannst Du mit einigermaßen kleinen Schritten einen Punkt in der Nähe suchen und dann mit einem der o.g. Billigverfahren recht schnell zu einer Lösung kommen.

Dann gibt es noch evolutionäre Programmierung/Optimierung, die immer funktioniert, aber eben manchmal ewig braucht.

Und andere numerische Verfahren, die ich nicht kenne, aber die bestimmt auch toll sind.

Schon mal gegoogelt?

[jfheins]

Du kannst die Methode der kleinsten Quadrate verwenden.
http://de.wikipedia.org/wiki/Methode...nsten_Quadrate
Damit kannst du eine nahezu beliebige Kurvengleichung an deine Punkte anpassen. Die Koeffizienten errechnen sich dann relativ einfach aus einem linearen Gleichungssystem.

Edit: Da das "relativ einfach" doch sehr relativ ist :stupid: habe ich dir mal eine Regression vorgerechnet.
Im Anhang sind die Maxima-Datei und ein Screenshot mit den Befehlen. Ich habe oben eine Beispielfunktion verwendet und ein bisschen an den Y-Werten gespielt und auf 2 Nachkommastellen angegeben. Am Ende kommen dann wieder ungefähr die Parameter 'raus.