math - systeme - système d'équation à 3 inconnues exercices corrigés




Résoudre une équation linéaire (7)

Êtes-vous à la recherche d'un progiciel permettant d'effectuer le travail ou effectuant réellement les opérations de la matrice, par exemple, et effectuant chaque étape?

Le premier, un de mes collègues vient d'utiliser Ocaml GLPK . C’est juste une enveloppe pour le GLPK , mais cela supprime beaucoup d’étapes de la configuration. Il semble que vous deviez vous en tenir au GLPK, en C, cependant. Pour ces derniers, grâce à Delicious pour avoir sauvegardé un vieil article, j’ai appris à faire de la LP, PDF . Si vous avez besoin d’aide spécifique pour la mise en place, faites-le nous savoir et, j'en suis sûr, nous reviendrons nous aider, mais je pense que la situation est assez simple. Bonne chance!

J'ai besoin de résoudre par programme un système d'équations linéaires en C, Objective C ou (si nécessaire) en C ++.

Voici un exemple des équations:

-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

A partir de cela, j'aimerais obtenir la meilleure approximation possible pour a , b et tx .


D'après le libellé de votre question, il semble que vous avez plus d'équations que d'inconnues et que vous souhaitez minimiser les incohérences. Cela se fait généralement avec une régression linéaire, qui minimise la somme des carrés des incohérences. En fonction de la taille des données, vous pouvez le faire dans un tableur ou dans un logiciel de statistiques. R est un package gratuit de haute qualité qui effectue une régression linéaire, parmi beaucoup d'autres choses. Il y a beaucoup de choses à la régression linéaire (et beaucoup de pièges), mais comme c'est simple à faire pour des cas simples. Voici un exemple R utilisant vos données. Notez que le "tx" est l'interception à votre modèle.

> 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  

Jetez un coup d'œil à Microsoft Solver Foundation .

Avec cela, vous pourriez écrire un code comme ceci:

  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); 

Voici la sortie:
=== Rapport de service Solver Foundation ===
Date / heure: 20/04/2009 23:29:55
Nom du modèle: par défaut
Capacités demandées: LP
Temps de résolution (ms): 1027
Durée totale (ms): 1414
Résoudre le statut d'achèvement: Optimal
Le solveur sélectionné: Microsoft.SolverFoundation.Solvers.SimplexSolver
Directives:
Microsoft.SolverFoundation.Services.Directive
Algorithme: Primal
Arithmétique: hybride
Prix ​​(exact): par défaut
Prix ​​(double): SteepestEdge
Base: mou
Nombre de pivots: 3
=== Détails de la solution ===
Buts:

Les décisions:
a: 0.0785250000000004
b: -0.180612500000001
c: -41.6375875


Personnellement, je suis partisan des algorithmes de Recettes numériques . (J'adore l'édition C ++.)

Ce livre vous apprendra pourquoi les algorithmes fonctionnent et vous montrera quelques implémentations plutôt bien déboguées de ces algorithmes.

Bien sûr, vous pouvez simplement utiliser aveuglément CLAPACK (je l’utilise avec beaucoup de succès), mais je voudrais d’abord taper à la main un algorithme d’élimination gaussien pour avoir au moins une petite idée du type de travail nécessaire à la création de ces algorithmes. stable.

Plus tard, si vous faites de l'algèbre linéaire plus intéressante, une recherche autour du code source d' Octave répondra à beaucoup de questions.


Vous pouvez résoudre ce problème avec un programme exactement de la même manière que vous le résolvez manuellement (avec multiplication et soustraction, puis en intégrant les résultats dans les équations). Ce sont de jolies mathématiques standard au niveau secondaire.

-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)

Donc vous vous retrouvez avec:

a =   0.0785250
b =  -0.1806125
c = -41.6375875

Si vous reconnectez ces valeurs à A, B et C, vous constaterez qu'elles sont correctes.

L'astuce consiste à utiliser une simple matrice 4x3 qui se réduit à une matrice 3x2, puis à un 2x1 qui est "a = n", n étant un nombre réel. Une fois que vous avez cela, vous l'insérez dans la matrice suivante pour obtenir une autre valeur, puis ces deux valeurs dans la matrice suivante jusqu'à ce que vous ayez résolu toutes les variables.

Si vous avez N équations distinctes, vous pouvez toujours résoudre pour N variables. Je dis distinct parce que ces deux ne sont pas:

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

Ce sont la même équation multipliée par deux, de sorte que vous ne pouvez pas obtenir de solution - multiplier la première par deux puis soustraire vous laisse avec la déclaration vraie mais inutile:

0 = 0 + 0

À titre d'exemple, voici un code C qui résout les équations simultanées que vous avez placées dans votre question. Tout d’abord quelques types nécessaires, des variables, une fonction de support pour l’impression d’une équation et le début de 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;

Ensuite, la réduction des trois équations à trois inconnues en deux équations à deux inconnues:

    // 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 ("");

Ensuite, la réduction des deux équations à deux inconnues en une équation à une inconnue:

    // 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 ("");

Maintenant que nous avons une formule du type number1 = unknown * number2 , nous pouvons simplement calculer la valeur unknown <- number1 / number2 avec unknown <- number1 / number2 . Ensuite, une fois que vous avez déterminé cette valeur, remplacez-la par l'une des équations à deux inconnues et calculez la seconde valeur. Puis substituez ces deux inconnus (maintenant connus) dans l'une des équations d'origine et vous avez maintenant les valeurs pour les trois inconnus:

    // 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;
}

La sortie de ce code correspond aux calculs précédents dans cette réponse:

         >: -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)

La règle de Cramer et l'élimination gaussienne sont deux bons algorithmes polyvalents (voir aussi Equations linéaires simultanées ). Si vous recherchez du code, consultez GiNaC , Maxima et SymbolicC++ (selon vos besoins en licences, bien sûr).

EDIT: Je sais que vous travaillez en C, mais je dois aussi parler de SymPy (un système de calcul SymPy en Python). Vous pouvez apprendre beaucoup de ses algorithmes (si vous pouvez lire un peu de python). En outre, il est sous la nouvelle licence BSD, alors que la plupart des packages mathématiques gratuits sont en GPL.


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






linear-equation