-
Notifications
You must be signed in to change notification settings - Fork 0
/
Matrix.pde
122 lines (100 loc) · 2.8 KB
/
Matrix.pde
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
class Matrix{
double[][] mat;
int rows, cols;
Matrix(int rows, int cols){
if(cols != rows+1 || rows <= 0) println("Invalid augmented matrix");
this.rows = rows;
this.cols = cols;
mat = new double[rows][cols];
}
void printe(){
for(int i=0; i<rows; i++){
for(int j=0; j<cols; j++){
print(mat[i][j]+" ");
}
print("\n");
}
}
// returns 1 if matrix is singular
int forwardElim(){
for (int k=0; k<rows; k++)
{
// init maximum value and index for pivot
int i_max = k;
double v_max = mat[i_max][k];
/* find greater amplitude for pivot if any */
for (int i = k+1; i < rows; i++){
if (Math.abs(mat[i][k]) > v_max){
v_max = mat[i][k];
i_max = i;
}
}
if (mat[k][i_max] == 0)
return 1; // Matrix is singular
// Swap the greatest value row with current row
if (i_max != k)
swap_row(k, i_max);
for (int i=k+1; i<rows; i++)
{
/* factor f to set current row kth elemnt to 0,
* and remaining kth column to 0 */
double f = mat[i][k]/mat[k][k];
/* subtract fth multiple of corresponding kth
row element*/
for (int j=k+1; j<cols; j++)
mat[i][j] -= mat[k][j]*f;
/* filling lower triangular matrix with zeros*/
mat[i][k] = 0;
}
//print();
}
//print();
return 0;
}
double[] backSub(){
double[] x = new double[rows]; // An array to store solution
/* Start calculating from last equation up to the
first */
for (int i = rows-1; i >= 0; i--)
{
/* start with the RHS of the equation */
x[i] = mat[i][rows];
/* Initialize j to i+1 since matrix is upper
triangular*/
for (int j=i+1; j<rows; j++)
{
/* subtract all the lhs values
* except the coefficient of the variable
* whose value is being calculated */
x[i] -= mat[i][j]*x[j];
}
/* divide the RHS by the coefficient of the
unknown being calculated */
x[i] = x[i]/mat[i][i];
}
double[] sol = new double[rows];
for (int i=0; i<rows; i++)
sol[i] = x[i];
return sol;
}
void swap_row(int i, int j){
for (int k=0; k<=rows; k++)
{
double temp = mat[i][k];
mat[i][k] = mat[j][k];
mat[j][k] = temp;
}
}
double[] gaussianElim(){
double[] solution;
int isSingular = this.forwardElim();
if(isSingular == 1){
println("Fucked up matrix, can't solve");
solution = new double[0];
}
else{
solution = backSub();
}
return solution;
}
}