Lineares Gleichungssystem lösen

[Kegasetu]

Hallo Zusammen.
Für mich und mein Programm geht es in die nächste Runde.

Ich stehe vor der Aufgabe Gleichungen aufzustellen und zu lösen.

Mein Überlegung ging in die Richtung Gauß. Da ich allerdings mir etwas Arbeit sparen wollte, will ich vorher die Frage loswerden, ab es dafür nicht eine vorgefertigte "solve"-Funktion gibt?

Im Internet habe ich nur Beiträge von vor 10 Jahren gefunden. Eventuell ist Delphi ja diesbezüglich etwas weiter.

[Rollo62]

Du könntest bei den Livebindings mal den Expression Parser ansehen,
aber ein fertiges Sover-Modul ist meines Wissens nicht dabei.
Bei TMS gibt es auch was, vermutlich aber auch nicht der fertige Solver.

[TigerLilly]

TMS kann das, die haben symbolische Mathe mit an Board - sehr cool. Aber zum Lösen linearer Gleichungssysteme mit Delphi/Pascal gibt Google auch einiges her.

Auch hier in der DP gibt es was:
https://www.delphipraxis.net/225-lin...en-loesen.html

[Kegasetu]

Ich werde mich mal reinlesen. Vielen Dank!

[Kegasetu]

Zitat:

Zitat von TigerLilly (Beitrag 1475899)

TMS kann das, die haben symbolische Mathe mit an Board - sehr cool. Aber zum Lösen linearer Gleichungssysteme mit Delphi/Pascal gibt Google auch einiges her.

Auch hier in der DP gibt es was:
https://www.delphipraxis.net/225-lin...en-loesen.html

Ich habe mich jetzt etwas mit TMS auseiandergesetzt, blicke aber nicht so ganz durch.. Info ist absolut nicht meine Stärke. Unter den Beispielen war einiges dabei, aber ich bin nicht so ganz durchgestiegen...

Gibt es zufällig ein Beispiel, bei dem ein solver in Delphi fertig implementiert ist? Natürlich nicht den ganzen Backgroud, sondern die Schreibweise aus der neuen Funktion heraus?

[Kegasetu]

Ich habe mich mit den Operatoren von TMS etwas länger auseinandergesetzt, jedoch liegt der Entschluss, dass dies nicht unbedingt notwendig ist. In Excel ist es relativ einfach solche Gleichungssysteme zu lösen, dafür gibt es fast schon komplett vorgefertigte Funktionen.
Während eine Internetsuche bin ich bezüglich Gauß und Delphi dann auf folgenden Code gestoßen:

Delphi-Quellcode:

Beispiel 15.2 Der allgemeine Gauß'sche Algorithmus.
Diesmal verwenden wir ein dynamische Array. Das zweidimensionale Array wird deklariert durch den Typ:

type TarrayOfArrayOfExtended = array of array of extended;

Falls koeff[i,i] = 0 ist müssen noch die Spalten vertauscht werden. Das geschieht mit der Permutation p und Umkehrpermutation q.
Hier musst Du noch die Untit MathMohr miteinbinden. (oder die Entsprechenden Prozeduren durch eigene ersetzten.)

function IsInteger(const x: extended; eps: extended): boolean;
begin //eps globale Variable. Zum Beispiel eps = 1E-9
  result := frac(abs(x) + eps) < eps*2;
end;

FUNCTION ggtInt(a,b:longint):longint;
  begin if b=0 then result:=a
        else result:=ggtInt(b,a mod b);
  end;

FUNCTION ggtReal(a,b:Extended):Extended;
  begin if (a < maxlongint) and (b < maxlongint) then
           result:=ggTInt(round(a),round(b))
        else Begin
          if abs(b) < 0.5 then result:=a
          else result:=ggtReal(b,a-b*int(a/b));
        end;
end;

FUNCTION kgVReal(a, b: extended): extended;
begin result := a * b / ggTReal(a, b) end;

procedure LSG(n: integer; var aa: TarrayOfArrayOfExtended; var xx: array of extended);
             //Arrray aa[0..n,0..n+1]
var i0, j0,j0max: integer; //n Unbekannte bzw. Gleichungen
  p, q: array of integer; //permutation q=inv_p
  procedure invers; //erzeugt q = inverse Permutation von p
  var u, v: integer;
  begin
    for u := 0 to n - 1 do
      for v := 0 to n - 1 do
        if p[u] = v then q[v] := u;
  end;
  procedure tausche_sp(i, j: integer); // Spalten werden ausgetauscht
  var u, k: integer;
    x: extended;
  begin
    for u := 0 to n - 1 do Begin
      x := aa[u, i];
      aa[u, i] := aa[u, j];
      aa[u, j] := x; //=altes aa[u,i]
    End;
    k := p[i];
    p[i] := p[j];
    p[j] := k; //altes p[i]
    invers;
  end;

  procedure macheZeileGanzzahlig(zeile: integer); //bis auf rechte Seite  aa[z,n+1]
  var k : integer;
    d, zae, nen: extended; //wird das kgV des Nenners
  begin
    try
      d := 1;
      for k := 0 to n - 1 do Begin
        if not ErmittleBruch(abs(aa[zeile, k]), zae, nen,g_eps) then exit;
        d := kgVReal(nen, d);
      End;
      for k := 0 to n do aa[zeile, k] := d * aa[zeile, k];
      //Jetzt noch kürzen
      if not isInteger(aa[zeile, 1], g_eps) then exit;
      d := round(aa[zeile, 1]);
      for k := 0 to n - 1 do Begin
        if not isInteger(aa[zeile, k], g_eps) then exit;
        if d = 0 then d := round(aa[zeile, k]); //falls aa[zeile,1..]=0
        if abs(aa[zeile, k]) > 0 then
          d := ggTReal(round(aa[zeile, k]), d);
      End;
      if d <> 0 then for k := 0 to n do aa[zeile, k] := aa[zeile, k] / d;
    except {dann halt nicht} end;
    if aa[zeile, zeile] < 0 then for k := 0 to n do aa[zeile, k] := -aa[zeile, k];
  end;

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

begin //Hauptprogramm
  setlength(p,n+1);
  setlength(q,n);
  for j0 := 0 to n - 1 do Begin
    p[j0] := j0;
    q[j0] := j0;
  End;
  for i0 := 0 to n - 1 do Begin
    j0max := i0;
    for j0 := i0 + 1 to n - 1 do if abs(aa[i0, j0]) > abs(aa[i0, j0max]) then j0max := j0;
    if aa[i0,j0max] = 0 then Begin
      showmessage('Keine eindeutige Lösng!');
      exit;
    End;
    VereinfacheRest(i0);
  end;
  for i0 := 0 to n - 1 do xx[p[i0]] := aa[i0, n] ;
end;

procedure ArrayInMemo(n: integer; aa: TarrayOfArrayOfExtended; m:Tmemo);
  var i, j: integer;   //aa n zeilen und n + 1 Spalten
      s: string;
begin
  m.Lines.Clear;
  for i := 0 to n-1 do Begin
    s := '';
    for j := 0 to n do
      s := s + ' ' + floatToStr(aa[i,j]);
    m.Lines.Add(s);
  End;
end;

procedure TForm1.Button1Click(Sender: TObject);
  var n, i: integer;
      koeff: TarrayOfArrayOfExtended;
      a: array of extended;
begin
  n := 3; // von 0 bis 2
  setlength(a, n);
  setlength(koeff, length(a),length(a) + 1);
  koeff[0,0] := 5;  koeff[0,1] := -7; koeff[0,2] := 1; koeff[0,3] := 9;
  koeff[1,0] := -3; koeff[1,1] := 2;  koeff[1,2] := 3; koeff[1,3] := 4;
  koeff[2,0] := 2;  koeff[2,1] := 3;  koeff[2,2] := 4; koeff[2,3] := 0;
  ArrayInMemo(length(a),Koeff,memo1);
  LSG(length(a),koeff,a);
  memo1.lines.add('Lösung');
  for i := 0 to n - 1 do
    memo1.lines.add(ReellZuBruch(a[i]));
end;

(https://kilchb.de/lektionen1ff.php)

Ich bekomme das Kotzen :D so kompliziert kann das doch nicht sein?

Gegeben ist eine 10x10 Matrix, welche gelöst werden muss.

[TigerLilly]

Bist du mit der Mathematik dahinter vertraut? Wenn nicht, ist das natürlich schwierig umzusetzen.

Aber sieh doch hier
https://www.delphipraxis.net/225-lin...en-loesen.html

Da ist auch die Mathe dahinter gut erklärt.

[Kegasetu]

Zitat:

Zitat von TigerLilly (Beitrag 1476341)

Bist du mit der Mathematik dahinter vertraut? Wenn nicht, ist das natürlich schwierig umzusetzen.

Aber sieh doch hier
https://www.delphipraxis.net/225-lin...en-loesen.html

Da ist auch die Mathe dahinter gut erklärt.

Die Mathematik dahinter ist ja nicht so kompliziert, aber beim Verständnis der Informatik hapert es gewaltig. Ich bin nur kurzfristig gezwungen mich damit auseinander zu setzen, also keine Sache die wirklich verinnerlichen muss.

Hier z.B.

Delphi-Quellcode:

type
  TGaussSolved = array of Extended;
  TGaussLine = TGaussSolved;
  TGaussMatrix = array of TGaussLine;

function SolveLinearSystem(A: TGaussMatrix; m, n: Integer): TGaussSolved;
var
  i, j, k: Integer;
  Pivot: TGaussLine;
  PivotRow: Integer;
  Multiplicator, Sum: Extended;

Was genau diese Type-Funktion macht, kann ich mir nicht ganz erklären. So etwas habe ich nie behandelt. Nicht falsch verstehen, ich bin nicht faul um mich darein zu denken, aber eine möglich einfachen Weg würde ich schon gerne gehen, da ich dies nur, wie schon erwähnt, einmal mache.

[Andreas13]

Hallo Kegasetu,
Die Parameter der Funktion sind wie folgt zu deuten:

Delphi-Quellcode:

A: TGaussMatrix;

ist die sog. erweiterte Koeffizienten-Matrix des Gleichungssystems, d.h. linke UND rechte Seiten aller Gleichungen zusammen:

a_0_0, a_0_1, a_0_2, ... a_0_m, b_0 // 1. Gleichung
a_1_0, a_1_1, a_1_2, ... a_1_m, b_1 // 2. Gleichung
... ... ... ... ... ... ... ... ... ... ... ... ... ... ...
a_m_0, a_m_1, a_m_2, ... a_m_m, b_m // letzte Gleichung

a_0_0 bis a_m_m sind die Koeffizienten der linken Seiten der Gleichungen.
Die Reihenfolge der Indizes der Matrix ist immer: Zeile, Spalte.
b_0 bis b_m sind die Koeffizienten der rechten Seiten der Gleichungen.
m: Anzahl der Zeilen (= Anzahl der Unbekannten)
n: Anzahl der Spalten
Es gilt hierbei: n = m + 1 d.h. n_Spalten = n_Zeilen + 1
Der Rückgabewert

Delphi-Quellcode:

TGaussSolved;

der Funktion ist der Lösungsvektor und enthält die berechneten Werte:
X_0, X_1, X_2, ... X_m

Die Vektoren und die Matrix sind hier jeweils NULL-basierte Arrays.
Viel Erfolg!
Gruß, Andreas

[Kegasetu]

Zitat:

Zitat von Andreas13 (Beitrag 1476361)

Hallo Kegasetu,
Die Parameter der Funktion sind wie folgt zu deuten:

Delphi-Quellcode:

A: TGaussMatrix;

ist die sog. erweiterte Koeffizienten-Matrix des Gleichungssystems, d.h. linke UND rechte Seiten aller Gleichungen zusammen:

a_0_0, a_0_1, a_0_2, ... a_0_m, b_0 // 1. Gleichung
a_1_0, a_1_1, a_1_2, ... a_1_m, b_1 // 2. Gleichung
... ... ... ... ... ... ... ... ... ... ... ... ... ... ...
a_m_0, a_m_1, a_m_2, ... a_m_m, b_m // letzte Gleichung

a_0_0 bis a_m_m sind die Koeffizienten der linken Seiten der Gleichungen.
Die Reihenfolge der Indizes der Matrix ist immer: Zeile, Spalte.
b_0 bis b_m sind die Koeffizienten der rechten Seiten der Gleichungen.
m: Anzahl der Zeilen (= Anzahl der Unbekannten)
n: Anzahl der Spalten
Es gilt hierbei: n = m + 1 d.h. n_Spalten = n_Zeilen + 1
Der Rückgabewert

Delphi-Quellcode:

TGaussSolved;

der Funktion ist der Lösungsvektor und enthält die berechneten Werte:
X_0, X_1, X_2, ... X_n

Die Vektoren und die Matrix sind hier jeweils NULL-basierte Arrays.
Viel Erfolg!
Gruß, Andreas

Danke, die Erklärung hilft.
Ich versuche das jetzt mal auf diesen Fall zu adaptieren:

l1*a+l2*b+l3*b+00*c+00*d+00*e+00*f+00*g+00*h+000*i =b1
00*a+l2*b+l3*b+l4*c+l5*d+00*e+00*f+00*g+00*h+000*i =b2
00*a+00*b+l3*b+l4*c+l5*d+l6*e+00*f+00*g+00*h+000*i =b3
00*a+00*b+00*b+l4*c+l5*d+l6*e+l7*f+00*g+00*h+000*i =b4
00*a+00*b+00*b+00*c+00*d+l6*e+l7*f+l8*g+00*h+000*i =b5
00*a+00*b+00*b+00*c+00*d+00*e+l7*f+l8*g+l9*h+000*i =b6
00*a+00*b+00*b+00*c+00*d+00*e+00*f+l8*g+l9*h+l10*i =b7
00*a+00*b+00*b+00*c+00*d+00*e+00*f+00*g+l9*h+l10*i =b8
00*a+00*b+00*b+00*c+00*d+00*e+00*f+00*g+00*h+l10*i =b9
00*a+00*b+00*b+00*c+00*d+00*e+00*f+00*g+00*h+000*i =b10