Lösen eines linearen Gleichungs

stimmen
35

Ich muß programmatisch ein System von linearen Gleichungen in C, Objective C lösen, oder (falls erforderlich) C ++.

Hier ist ein Beispiel der Gleichungen:

-44.3940 = a * 50.0 + b * 37.0 + tx
-45.3049 = a * 43.0 + b * 39.0 + tx
-44.9594 = a * 52.0 + b * 41.0 + tx

Daraus würde Ich mag für die beste Annäherung zu erhalten a, bund tx.

Veröffentlicht am 03/08/2008 um 19:14
quelle vom benutzer
In anderen Sprachen...                            


10 antworten

stimmen
17

Cramer-Regel und Gauß - Elimination sind zwei gute, universelle Algorithmen (siehe auch lineare Gleichungssysteme ). Wenn Sie Code suchen, überprüfen GiNaC , Maxima und SymbolicC ++ (je nach Lizenzanforderungen, natürlich).

EDIT: Ich weiß , dass Sie in C Land arbeiten, aber ich muss auch für ein gutes Wort setzen SymPy (ein Computer - Algebra - System in Python). Sie können eine Menge von seinen Algorithmen lernen (wenn Sie ein wenig Python lesen kann). Außerdem ist es unter dem neuen BSD - Lizenz, während die meisten der kostenlosen Mathe - Pakete GPL sind.

Beantwortet am 03/08/2008 um 19:37
quelle vom benutzer

stimmen
15

Sie können dies mit einem Programm lösen genau die gleiche Art und Weise Sie es von Hand lösen (mit Multiplikation und Subtraktion, dann Ergebnisse zurück in die Gleichungen Fütterung). Das ist ziemlich Standard Gymnasialebene Mathematik.

-44.3940 = 50a + 37b + c (A)
-45.3049 = 43a + 39b + c (B)
-44.9594 = 52a + 41b + c (C)

(A-B): 0.9109 =  7a -  2b (D)
(B-C): 0.3455 = -9a -  2b (E)

(D-E): 1.2564 = 16a (F)

(F/16):  a = 0.078525 (G)

Feed G into D:
       0.9109 = 7a - 2b
    => 0.9109 = 0.549675 - 2b (substitute a)
    => 0.361225 = -2b (subtract 0.549675 from both sides)
    => -0.1806125 = b (divide both sides by -2) (H)

Feed H/G into A:
       -44.3940 = 50a + 37b + c
    => -44.3940 = 3.92625 - 6.6826625 + c (substitute a/b)
    => -41.6375875 = c (subtract 3.92625 - 6.6826625 from both sides)

So können Sie am Ende mit:

a =   0.0785250
b =  -0.1806125
c = -41.6375875

Wenn Sie diese Werte wieder in A, B und C stecken, finden Sie sie richtig sind.

Der Trick besteht darin, eine einfache 4x3-Matrix zu verwenden, die auf eine 3x2-Matrix wiederum verringert, wird eine 2x1, die „a = n“, n eine tatsächliche Zahl ist. Sobald Sie, dass Sie es in die nächste Matrix aufzupäppeln einen anderen Wert zu erhalten, dann diese beiden Werte in die oben nächste Matrix, bis Sie alle Variablen gelöst haben.

Vorausgesetzt, man N verschiedene Gleichungen haben, können Sie immer für N Variablen lösen. Ich sage anders, weil diese beiden sind nicht:

 7a + 2b =  50
14a + 4b = 100

Sie sind die gleiche Gleichung mit zwei multipliziert , so dass Sie sich nicht um eine Lösung von ihnen erhalten - Multiplizieren der ersten von zwei dann Subtrahieren Blätter Sie mit der wahren , aber nutzloser Aussage:

0 = 0 + 0

Als Beispiel ist hier einig C - Code, der die simultanen Gleichungen ausarbeitet , die Sie in Ihrer Frage gestellt sind. Zunächst einige notwendige Typen, Variablen, eine Stützfunktion für eine Gleichung auszudrucken, und der Anfang main:

#include <stdio.h>

typedef struct { double r, a, b, c; } tEquation;
tEquation equ1[] = {
    { -44.3940,  50, 37, 1 },      // -44.3940 = 50a + 37b + c (A)
    { -45.3049,  43, 39, 1 },      // -45.3049 = 43a + 39b + c (B)
    { -44.9594,  52, 41, 1 },      // -44.9594 = 52a + 41b + c (C)
};
tEquation equ2[2], equ3[1];

static void dumpEqu (char *desc, tEquation *e, char *post) {
    printf ("%10s: %12.8lf = %12.8lfa + %12.8lfb + %12.8lfc (%s)\n",
        desc, e->r, e->a, e->b, e->c, post);
}

int main (void) {
    double a, b, c;

Als nächstes wird die Reduktion der drei Gleichungen mit drei Unbekannten zu zwei Gleichungen mit zwei Unbekannten:

    // First step, populate equ2 based on removing c from equ.

    dumpEqu (">", &(equ1[0]), "A");
    dumpEqu (">", &(equ1[1]), "B");
    dumpEqu (">", &(equ1[2]), "C");
    puts ("");

    // A - B
    equ2[0].r = equ1[0].r * equ1[1].c - equ1[1].r * equ1[0].c;
    equ2[0].a = equ1[0].a * equ1[1].c - equ1[1].a * equ1[0].c;
    equ2[0].b = equ1[0].b * equ1[1].c - equ1[1].b * equ1[0].c;
    equ2[0].c = 0;

    // B - C
    equ2[1].r = equ1[1].r * equ1[2].c - equ1[2].r * equ1[1].c;
    equ2[1].a = equ1[1].a * equ1[2].c - equ1[2].a * equ1[1].c;
    equ2[1].b = equ1[1].b * equ1[2].c - equ1[2].b * equ1[1].c;
    equ2[1].c = 0;

    dumpEqu ("A-B", &(equ2[0]), "D");
    dumpEqu ("B-C", &(equ2[1]), "E");
    puts ("");

Als nächstes wird die Reduktion der beiden Gleichungen mit zwei Unbekannten zu einer Gleichung mit einer Unbekannten:

    // Next step, populate equ3 based on removing b from equ2.

    // D - E
    equ3[0].r = equ2[0].r * equ2[1].b - equ2[1].r * equ2[0].b;
    equ3[0].a = equ2[0].a * equ2[1].b - equ2[1].a * equ2[0].b;
    equ3[0].b = 0;
    equ3[0].c = 0;

    dumpEqu ("D-E", &(equ3[0]), "F");
    puts ("");

Nun , da wir eine Formel des Typs haben number1 = unknown * number2, können wir einfach ausrechnen den unbekannten Wert mit unknown <- number1 / number2. Dann, sobald Sie aus diesem Wert rechneten haben, ersetzen Sie es in einer der Gleichungen mit zwei Unbekannten und arbeitet , um den zweiten Wert aus. Dann ersetzen diese beiden (jetzt bekannt) Unbekannten in einer der ursprünglichen Gleichungen und Sie nun die Werte für alle drei Unbekannten haben:

    // Finally, substitute values back into equations.

    a = equ3[0].r / equ3[0].a;
    printf ("From (F    ), a = %12.8lf (G)\n", a);

    b = (equ2[0].r - equ2[0].a * a) / equ2[0].b;
    printf ("From (D,G  ), b = %12.8lf (H)\n", b);

    c = (equ1[0].r - equ1[0].a * a - equ1[0].b * b) / equ1[0].c;
    printf ("From (A,G,H), c = %12.8lf (I)\n", c);

    return 0;
}

Der Ausgang dieses Codes entspricht die früheren Berechnungen in dieser Antwort:

         >: -44.39400000 =  50.00000000a +  37.00000000b +   1.00000000c (A)
         >: -45.30490000 =  43.00000000a +  39.00000000b +   1.00000000c (B)
         >: -44.95940000 =  52.00000000a +  41.00000000b +   1.00000000c (C)

       A-B:   0.91090000 =   7.00000000a +  -2.00000000b +   0.00000000c (D)
       B-C:  -0.34550000 =  -9.00000000a +  -2.00000000b +   0.00000000c (E)

       D-E:  -2.51280000 = -32.00000000a +   0.00000000b +   0.00000000c (F)

From (F    ), a =   0.07852500 (G)
From (D,G  ), b =  -0.18061250 (H)
From (A,G,H), c = -41.63758750 (I)
Beantwortet am 26/02/2009 um 11:59
quelle vom benutzer

stimmen
7

Für ein 3x3-System von linearen Gleichungen ich denke, es wäre in Ordnung Ihre eigenen Algorithmen ausrollen.

Allerdings könnten Sie Genauigkeit, Division durch Null oder wirklich kleine Zahlen kümmern und was ist unendlich viele Lösungen zu tun. Mein Vorschlag ist , mit einem Standard - numerischen linearen Algebra - Paket wie gehen LAPACK .

Beantwortet am 08/08/2008 um 18:59
quelle vom benutzer

stimmen
6

Werfen Sie einen Blick auf die Microsoft Solver Foundation .

Mit ihm können Sie Code wie folgt schreiben:

  SolverContext context = SolverContext.GetContext();
  Model model = context.CreateModel();

  Decision a = new Decision(Domain.Real, "a");
  Decision b = new Decision(Domain.Real, "b");
  Decision c = new Decision(Domain.Real, "c");
  model.AddDecisions(a,b,c);
  model.AddConstraint("eqA", -44.3940 == 50*a + 37*b + c);
  model.AddConstraint("eqB", -45.3049 == 43*a + 39*b + c);
  model.AddConstraint("eqC", -44.9594 == 52*a + 41*b + c);
  Solution solution = context.Solve();
  string results = solution.GetReport().ToString();
  Console.WriteLine(results); 

Hier ist die Ausgabe:
=== Solver Foundation Service Report === für
Datum und Uhrzeit: 04/20/2009 23.29.55
Modellname:
Standardfunktionen angefordert: LP
Zeit lösen (ms): 1027
Gesamtspielzeit (ms): 1414
Solve Completion Status: Optimal
Solver ausgewählt: Microsoft.SolverFoundation.Solvers.SimplexSolver
Richtlinien:
Microsoft.SolverFoundation.Services.Directive
Algorithmus: Primal
Arithmetik: Hybrid
Pricing (exakt): Standard
Pricing (double): SteepestEdge
Basis: Slack
Pivot Anzahl: 3
== = Lösungsdetails ===
Ziele:

Entscheidungen:
a: ,0785250000000004
b: -,180612500000001
c: -41,6375875

Beantwortet am 21/04/2009 um 04:33
quelle vom benutzer

stimmen
3

Template Numerical Toolkit von NIST hat Tools für das zu tun.

Einer der zuverlässigere Möglichkeiten ist die Verwendung eines QR - Zerlegung .

Hier ist ein Beispiel für einen Wrapper, so dass ich anrufen kann „GetInverse (A InvG)“ in meinem Code und es wird die inverse in InvG setzen.

void GetInverse(const Array2D<double>& A, Array2D<double>& invA)
   {
   QR<double> qr(A);  
   invA = qr.solve(I); 
   }

Array2D wird in der Bibliothek definiert.

Beantwortet am 25/08/2008 um 19:53
quelle vom benutzer

stimmen
3

Suchen Sie ein Softwarepaket, das die Arbeit oder tatsächlich tun, die Matrix-Operationen und so und machen jeden Schritt tun würde?

Die die erste, ein Kollege von mir verwendet nur Ocaml GLPK . Es ist nur ein Wrapper für die GLPK , aber es entfernt viele der Schritte aufzusetzen. Es sieht aus wie Sie gehen mit der GLPK zu haften zu haben, in C, though. Für letztere dank köstlich für einen alten Artikel Speicher verwendeten I LP eine Weile zurück, lernen PDF . Wenn Sie spezifische Hilfe benötigen weitere Einrichtung, lassen Sie uns wissen , und ich bin sicher, ich oder jemand in wandern zurück und helfen, aber ich denke , es ist ziemlich einfach von hier. Viel Glück!

Beantwortet am 03/08/2008 um 19:27
quelle vom benutzer

stimmen
2

In Bezug auf die Laufzeiteffizienz, andere haben beantwortet besser als ich Wenn Sie immer die gleiche Anzahl von Gleichungen haben als Variablen, Ich mag Cramer-Regel , da es einfach zu implementieren ist. Schreiben Sie einfach eine Funktion Determinante einer Matrix zu berechnen (oder eine , die schon geschrieben, ich bin sicher , dass Sie ein da draußen finden), und teilt die Determinanten zweier Matrizen.

Beantwortet am 19/09/2008 um 04:36
quelle vom benutzer

stimmen
2

Aus dem Wortlaut der Frage, so scheint es, wie Sie mehr Gleichungen als Unbekannte haben und Sie die Inkonsistenzen zu minimieren. Dies wird in der Regel mit linearer Regression durchgeführt, die die Summe der Quadrate der Abweichungen minimieren. Je nach Größe der Daten, können Sie dies in einer Tabelle tun oder in einem statistischen Paket. R ist ein hochwertiges, freies Paket, das die lineare Regression der Fall ist, unter vielen anderen Dingen. Es gibt eine Menge zu linearer Regression (und vielen Gotcha Jahre), aber es ist einfach für einfache Fälle zu tun. Hier ist ein Beispiel R Ihrer Daten. Beachten Sie, dass die „TX“ ist der Schnittpunkt zu Ihrem Modell.

> y <- c(-44.394, -45.3049, -44.9594)
> a <- c(50.0, 43.0, 52.0)
> b <- c(37.0, 39.0, 41.0)
> regression = lm(y ~ a + b)
> regression

Call:
lm(formula = y ~ a + b)

Coefficients:
(Intercept)            a            b  
  -41.63759      0.07852     -0.18061  
Beantwortet am 17/09/2008 um 15:22
quelle vom benutzer

stimmen
1
function x = LinSolve(A,y)
%
% Recursive Solution of Linear System Ax=y
% matlab equivalent: x = A\y 
% x = n x 1
% A = n x n
% y = n x 1
% Uses stack space extensively. Not efficient.
% C allows recursion, so convert it into C. 
% ----------------------------------------------
n=length(y);
x=zeros(n,1);
if(n>1)
    x(1:n-1,1) = LinSolve( A(1:n-1,1:n-1) - (A(1:n-1,n)*A(n,1:n-1))./A(n,n) , ...
                           y(1:n-1,1) - A(1:n-1,n).*(y(n,1)/A(n,n))); 
    x(n,1) = (y(n,1) - A(n,1:n-1)*x(1:n-1,1))./A(n,n); 
else
    x = y(1,1) / A(1,1);
end
Beantwortet am 30/12/2014 um 15:24
quelle vom benutzer

stimmen
1

Ich persönlich bin teilweise auf die Algorithmen der Numerical Recipes . (Ich bin begeistert von der C ++ Ausgabe.)

Dieses Buch wird Ihnen zeigen, warum die Algorithmen arbeiten, und zeigen Ihnen einige ziemlich gut-ausgetestet Implementierungen dieser Algorithmen.

Natürlich könnte man einfach blind verwenden CLAPACK (ich es mit großem Erfolg verwendet haben), aber ich würde eine schwache Vorstellung von der Art der Arbeit eine Gaußsche Eliminationsalgorithmus zumindest erste Hand geben, die in der Herstellung dieser Algorithmen gegangen ist stabil.

Später, wenn Sie mehr interessante lineare Algebra tun, um den Quellcode sucht Octave wird eine Menge Fragen beantworten.

Beantwortet am 25/08/2008 um 19:22
quelle vom benutzer

Cookies help us deliver our services. By using our services, you agree to our use of cookies. Learn more