[UVA][高斯消去][線性系統] 10109 - Solving Systems of Linear ...
Solving Systems of Linear Equations |
The Problem
You may have solved linear equation early in the school. Problems involving solving sets of linear equation are very important in the field of Engineering and Mathematics.Let us consider that we have a system of linear equations
a11x1 + a12x2 + a13x3 = c1
a21x1 + a22x2 + a23x3 = c2
a31x1 + a32x2 + a33x3 = c3
We can solve it by reducing technique:
Step 1:
x1 + a12/a11x2 + a13/a11x3 = c1/a11
a21x1 + a22x2 + a23x3 = c2
a31x1 + a32x2 + a33x3 = c3
Step
2:
x1 + a12/a11x2 + a13/a11x3
= c1/a11
(-a21* a12/a11)a22x2
+ (-a21* a13/a11)a23x3
= (-a21* c1/a11)c2
(-a31* a12/a11)a22x2
+ (-a31* a13/a11)a23x3
= (-a31* c1/a11)c2
Now do as step 1 for second row and so on.
This can be made more effective using matrix method. The set of equation
for `n' unknowns can be written as
a11x1 + a12x2 + a13x3
+ ... + a1nxn = c1
a21x1 + a22x2 + a23x3
+ ... + a2nxn = c2
a31x1 + a32x2 + a33x3
+ ... + a3nxn = c3
... ... ... ... ...
... ... ... ... ...
... ... ... ...
am1x1 + am2x2 + am3x3
+ ... + amnxn = cn
In matrix form
Compactly [A] * {X} = {C}
From this we can solve values of X's. The matrix [AC] is called an
augmented (see example below) matrix. If after elimination process the
rank of matrix [A] and rank of matrix [AC] not equals, the system is called
inconsistent and it does not have a solution. If the matrix is consistent
and number of unknowns is greater then rank of matrix then the matrix system
has arbitarily many solutions containing (NumberOfUnknowns-rank) arbitary
constants. Rank of a matrix is defined as the number of non zero rows of
a matrix system. Otherwise if the rank and number of unknows equals then
the system has been solved.
For example let a system of equations be
9x1 + 4x2 + x3 = -17
x1 - 2x2 - 6x3 = 14
x1 + 6x2 = 4
This sets of equation can be written as
9 4 1 -17
1 -2 -6 14
1 6 0
4
So the steps involving the solution is
Step : 1
1
-2
-6
14
1
6
0
4
9
4
1 -17
Step : 2
1
-2
-6
14
0
8
6 -10
0
22
55 -143
Step : 3
1
-2
-6
14
0
1
3/4
-5/4
0
0
77/2 -231/2
Step : 4
1
-2
-6
14
0
1
3/4
-5/4
0
0
1
-3
Step : 5
1
-2
0
-4
0
1
0
1
0
0
1
-3
Step : 6
1
0
0
-2
0
1
0
1
0
0
1
-3
x1 = -2
x2 = 1
x3 = -3
Again consider this system
2x1 + 2x2
+ 2x3 = 2
4x1 + 4x2
+ 4x3 = 4
16x1 + 16x2
+ 16x3 = 16
Steps are:
Step : 1
2
2
2
2
4
4
4
4
16
16
16
16
Step : 2
1
1
1
1
0
0
0
0
0
0
0
0
This system has number of unknowns 3 and rank is 1. So this system
has arbitarily many solutions containing (3-1) = 2 arbitary constants.
Another system
x + y = 10
x + y = 20
x + y = 50
Steps are:
Step : 1
1
1
10
1
1
20
2
2
50
Step : 2
1
1
10
0
0
10
0
0
30
Step : 3
1
1
10
0
0
10
0
0
30
Step : 4
1
1
10
0
0
10
0
0
1
Step : 5
1
1
0
0
0
0
0
0
1
Step : 6
1
1
0
0
0
0
0
0
1
As rank of [A] (in this case: rank(A) = 1) is not equal to the rank of augmented matrix [AC] (in this case: rank(AC) = 2) , the system has no solution.
However though there are other methods to compute this solution for
the matrix system, the main problem occurs are
1. Round off errors or computational
error due to the use of floating point number
2. Error due to wrong order
of the given equation.
To prevent round off error due to floating point number an approach can be used, similar to the process of doing fractional number. So we may use 1/3 as a expression of two integer, the numerator and the denominator, instead of .333333 (with loss of precision). Thus we can prevent this kind of error.
Consider this set of equations
5x3 = 10
3x2 - 3x3 = 3
2x1 - x2 + 2x3 = 7
This set of equations can be written as
0 0 5
10
0 3 -3
3
2 -1 2
7
Now how will you evaluate this matrix without ordering?
The Input
The first line of input is the number of the problem. The next line contains two integers - NumberOfUnknowns and NumberOfEquations (none of these is less then or equal to 0 and greater then 50). The next lines contains the matrix for the system of linear equations. There are number of rows equal to the NumberOfEquations and number of column equal to the NumberOfUnknowns+1. The numbers may be fractional, that is there may be numbers like 1/3 or 6/8. An problem number zero indicates the end of input.The Output
First print (without the quotation mark)"Solution for Matrix System # N"Here `N' is the problem number as taken from input. Then on the next line, for each system of equations output the solution (if exists) expressed in the fractional form in each line. You may assume each of the numerator and denominator part will not exceed the limit of data type long long (64 bit). If there are many solutions as described above print (without the quotation mark)
"Infinitely many solutions containing n arbitrary constants."(here `n' is the number as described above) , and if there is no solutions print (without the quotation mark)
"No Solution.".Print a blank line between two systems of linear equations.
Sample Input
1 3 3 9 4 1 -17 1 -2 -6 14 1 6 0 4 2 3 3 2 2 2 2 4 4 4 4 16/1 16/1 16/1 16/1 3 2 3 1 1 10 1 1 20 2 2 50 4 1 1 3 10 0
Sample Output
Solution for Matrix System # 1 x[1] = -2 x[2] = 1 x[3] = -3 Solution for Matrix System # 2 Infinitely many solutions containing 2 arbitrary constants. Solution for Matrix System # 3 No Solution. Solution for Matrix System # 4 x[1] = 10/3
Suman Mahbub
01-04-2001
#include <cstdio>
#include <iostream>
#include <algorithm>
using namespace std;
class Frac {
public:
long long a, b;
Frac() {
a = 0, b = 1;
}
Frac(int x, int y) {
a = x, b = y;
reduce();
}
Frac operator+(const Frac &y) {
long long ta, tb;
tb = this->b/gcd(this->b, y.b)*y.b;
ta = this->a*(tb/this->b) + y.a*(tb/y.b);
Frac z(ta, tb);
return z;
}
Frac operator-(const Frac &y) {
long long ta, tb;
tb = this->b/gcd(this->b, y.b)*y.b;
ta = this->a*(tb/this->b) - y.a*(tb/y.b);
Frac z(ta, tb);
return z;
}
Frac operator*(const Frac &y) {
long long tx, ty, tz, tw, g;
tx = this->a, ty = y.b;
g = gcd(tx, ty), tx /= g, ty /= g;
tz = this->b, tw = y.a;
g = gcd(tz, tw), tz /= g, tw /= g;
Frac z(tx*tw, ty*tz);
return z;
}
Frac operator/(const Frac &y) {
long long tx, ty, tz, tw, g;
tx = this->a, ty = y.a;
g = gcd(tx, ty), tx /= g, ty /= g;
tz = this->b, tw = y.b;
g = gcd(tz, tw), tz /= g, tw /= g;
Frac z(tx*tw, ty*tz);
return z;
}
private:
static long long gcd(long long x, long long y) {
if(!y) return x;
if(x < 0) x *= -1;
if(y < 0) y *= -1;
long long t;
while(x%y)
t = x, x = y, y = t%y;
return y;
}
void reduce() {
long long g = gcd(a, b);
a /= g, b /= g;
if(b < 0) a *= -1, b *= -1;
}
};
ostream& operator<<(ostream& out, const Frac&x) {
out << x.a;
if(x.b != 1)
out << '/' << x.b;
return out;
}
int main() {
int n, m, i, j, k, N;
char NUM[100], first = 0;
long long X, Y;
Frac matrix[100][100];
while(scanf("%d", &N) == 1 && N) {
scanf("%d %d", &n, &m);
for(i = 0; i < m; i++) {
for(j = 0; j <= n; j++) {
scanf("%s", NUM);
if(sscanf(NUM, "%lld/%lld", &X, &Y) == 2) {
matrix[i][j].a = X;
matrix[i][j].b = Y;
} else
sscanf(NUM, "%lld", &matrix[i][j].a), matrix[i][j].b = 1;
}
}
Frac tmp, one(1,1);
int idx = 0, rank = 0;
for(i = 0; i < m; i++) {
while(idx < n) {
int ch = -1;
for(j = i; j < m; j++)
if(matrix[j][idx].a) {
ch = j;
break;
}
if(ch == -1) {
idx++;
continue;
}
if(i != ch)
for(j = idx; j <= n; j++)
swap(matrix[ch][j], matrix[i][j]);
break;
}
if(idx >= n) break;
tmp = one/matrix[i][idx];
rank++;
for(j = idx; j <= n; j++)
matrix[i][j] = matrix[i][j]*tmp;
for(j = 0; j < m; j++) {
if(i == j) continue;
tmp = matrix[j][i]/matrix[i][idx];
for(k = idx; k <= n; k++) {
matrix[j][k] = matrix[j][k] - tmp*matrix[i][k];
}
}
idx++;
}
if(first) puts("");
first = 1;
printf("Solution for Matrix System # %d\n", N);
int sol = 0;
for(i = 0; i < m; i++) {
for(j = 0; j < n; j++) {
if(matrix[i][j].a)
break;
}
if(j == n) {
if(matrix[i][n].a == 0 && sol != 1)
sol = 0; // INFINITELY
else
sol = 1; // No Solution.
}
}
if(rank == n && sol == 0) {
for(i = 0; i < n; i++) {
printf("x[%d] = ", i+1);
cout << matrix[i][n] << endl;
}
continue;
}
if(sol == 1)
puts("No Solution.");
else
printf("Infinitely many solutions containing %d arbitrary constants.\n", n-rank);
}
return 0;
}