Member Avatar for 9tontruck

Hi..
d1 = a1x1 + b1x2 + c1x3
d2 = a2x1 + b2x2 + c2x3
d3 = a3x1 + b3x2 + c3*x3

Knowing all of a,b,c,d, I have to find x1, x2, x3. As you might noticed, this is high school math. But, how do I write the code to solve this?

/*
/a1 b1 c1| /x1| /d1|
|a2 b2 c2|*|x2|=|d2|
|a3 b3 c3/ |x3/ |d3|
*/
double A[3][3], X[3], D[3];
X[0] = ?
X[1] = ?
X[2] = ?

Can anyone give me an idea please?

Edited by pritaeas because: Moved.
  • 2 Contributors
  • 2 Replies
  • 426 Views
  • 5 Hours Discussion Span
  • Latest Post Latest Post by DavidB

Dani AI

Generated

Short answer: do Gaussian elimination with partial pivoting and back-substitution rather than explicitly forming the inverse. is right that triangularising the matrix and back-substituting is the standard, stable approach; and while showed the inverse-by-determinant route, that is more error prone in floating-point and unnecessary for a 3x3 solver.

The idea in words: form the 3x3 augmented matrix [A|d], for k = 0..2 pick the row with largest absolute value in column k and swap it into row k (partial pivoting). If the pivot is essentially zero (|pivot| < eps) the system is singular or ill-conditioned. Otherwise eliminate below the pivot to make a triangular system, then compute x by back-substitution. Check the residual afterwards to catch numerical issues.

A compact C++ implementation (self-contained for 3x3):

#include <cmath>

bool solve3x3(double A[3][3], double b[3], double x[3]) {
    double aug[3][4];
    for (int i=0;i<3;i++) {
        for (int j=0;j<3;j++) aug[i][j] = A[i][j];
        aug[i][3] = b[i];
    }
    const double eps = 1e-12;
    for (int k=0;k<3;k++) {
        int piv = k;
        double maxv = fabs(aug[k][k]);
        for (int i=k+1;i<3;i++) {
            double v = fabs(aug[i][k]);
            if (v > maxv) { maxv = v; piv = i; }
        }
        if (piv != k) for (int j=k;j<4;j++) std::swap(aug[k][j], aug[piv][j]);
        if (fabs(aug[k][k]) < eps) return false; // singular / no unique solution
        for (int i=k+1;i<3;i++) {
            double f = aug[i][k] / aug[k][k];
            for (int j=k;j<4;j++) aug[i][j] -= f * aug[k][j];
        }
    }
    for (int i=2;i>=0;i--) {
        double s = aug[i][3];
        for (int j=i+1;j<3;j++) s -= aug[i][j] * x[j];
        x[i] = s / aug[i][i];
    }
    return true;
}

After solving, compute the residual r = A*x - d and check its max absolute component (tolerance ~1e-9). For production or larger systems use a tested library (Eigen, GSL) to avoid hand-rolled pitfalls.

Recommended Answers

All 2 Replies

Member Avatar for 9tontruck 
public static function solveMatrix(a:Array, b:Array, X:Array):void
        {
            var invA:Array = [];
            invA[0] = [];
            invA[1] = [];
            invA[2] = [];

             var invDet:Number = 1/((a[0][0]*(a[2][2]*a[1][1]-a[2][1]*a[1][2]))-(a[1][0]*(a[2][2]*a[0][1]-a[2][1]*a[0][2]))+(a[2][0]*(a[1][2]*a[0][1]-a[1][1]*a[0][2])));
            invA[0][0] = invDet * (a[2][2]*a[1][1]-a[2][1]*a[1][2]);
            invA[0][1] = invDet * (-1 * (a[2][2]*a[0][1]-a[2][1]*a[0][2]));
            invA[0][2] = invDet * (a[1][2]*a[0][1]-a[1][1]*a[0][2]);

            invA[1][0] = invDet * (-1 * (a[2][2]*a[1][0]-a[2][0]*a[1][2]));
            invA[1][1] = invDet * (a[2][2]*a[0][0]-a[2][0]*a[0][2]);
            invA[1][2] = invDet * (-1 * (a[1][2]*a[0][0]-a[1][0]*a[0][2]));

            invA[2][0] = invDet * (a[2][1]*a[1][0]-a[2][0]*a[1][1]);
            invA[2][1] = invDet * (-1 * (a[2][1]*a[0][0]-a[2][0]*a[0][1]));
            invA[2][2] = invDet * (a[1][1]*a[0][0]-a[1][0]*a[0][1]);

            X[0] = invA[0][0] * b[0] + invA[0][1] * b[1] + invA[0][2] * b[2];
            X[1] = invA[1][0] * b[0] + invA[1][1] * b[1] + invA[1][2] * b[2];
            X[2] = invA[2][0] * b[0] + invA[2][1] * b[1] + invA[2][2] * b[2];
        }
Member Avatar for DavidB

You are doing this in high school?! It must be a pretty advanced high school. I don't recall doing anything like this before first year university.

In any case, normally, to solve N Equations in N Unknowns, [A](x) = (b), the easiest way is to make the [A] matrix triangular, and then use back-substitution. However, I doubt you want anything that complicated.

Are you only going to be working with systems of size 3? Are you supposed to solve the system using determinants? It looks like you are off to a good start: computing the determinant by the standard formula.

Are you unclear about the algorithm? Or the coding?

Edited by DavidB
Be a part of the DaniWeb community

We're a friendly, industry-focused community of developers, IT pros, digital marketers, and technology enthusiasts meeting, networking, learning, and sharing knowledge.