71 lines
2.1 KiB
C++
71 lines
2.1 KiB
C++
#include "Solver.h"
|
|
|
|
#include "Gauss.h"
|
|
|
|
static int FirstNotNullElementIndexOnLine(const Matrix& mat, std::size_t line) {
|
|
for (std::size_t i = 0; i < mat.GetColumnCount(); i++) {
|
|
if (!IsEqualZero(mat.at(line, i))) {
|
|
return i;
|
|
}
|
|
}
|
|
return -1;
|
|
}
|
|
|
|
Vect Solver::Image(Matrix&& a_Matrix) const {
|
|
a_Matrix.Transpose();
|
|
Gauss::GaussJordan(a_Matrix, false, false);
|
|
a_Matrix.Transpose();
|
|
return {std::move(a_Matrix)};
|
|
}
|
|
|
|
// https://en.wikipedia.org/wiki/Kernel_(linear_algebra)#Computation_by_Gaussian_elimination
|
|
Vect Solver::Kernel(Matrix&& a_Matrix) const {
|
|
std::size_t matrixRawCount = a_Matrix.GetRawCount();
|
|
std::size_t matrixColumnCount = a_Matrix.GetColumnCount();
|
|
|
|
a_Matrix.Transpose();
|
|
a_Matrix.Augment(Matrix::Identity(a_Matrix.GetRawCount()));
|
|
Gauss::GaussJordan(a_Matrix, false, true);
|
|
a_Matrix.Transpose();
|
|
|
|
// nombre de colonnes non nulles
|
|
std::size_t origine_colonne = Vect(a_Matrix.SubMatrix(0, 0, matrixRawCount, matrixColumnCount)).GetCardinal();
|
|
|
|
return {a_Matrix.SubMatrix(
|
|
matrixRawCount, origine_colonne, a_Matrix.GetRawCount() - matrixRawCount, a_Matrix.GetColumnCount() - origine_colonne)};
|
|
}
|
|
|
|
VectAffine Solver::RectangularSystem(Matrix&& a_MatrixA, const Matrix& a_VectorB) const {
|
|
Matrix mat = a_MatrixA;
|
|
mat.Augment(a_VectorB);
|
|
Gauss::GaussJordan(mat, true, true);
|
|
|
|
Solver solver;
|
|
|
|
Vect noyau = solver.Kernel(std::move(a_MatrixA));
|
|
Matrix origin = mat.SubMatrix(0, mat.GetColumnCount() - 1, mat.GetRawCount(), 1);
|
|
|
|
// on calcule le vecteur qui dirige l'espace affine
|
|
Matrix fullOrigin {mat.GetColumnCount() - 1, 1};
|
|
for (std::size_t i = 0; i < mat.GetRawCount(); i++) {
|
|
int pivot_index = FirstNotNullElementIndexOnLine(mat, i);
|
|
|
|
if (static_cast<std::size_t>(pivot_index) == mat.GetColumnCount() - 1) {
|
|
// on a une ligne du type 0 = n. Aucune solution !
|
|
return {Matrix {}, Matrix::ColumnVector({0})};
|
|
}
|
|
|
|
// ligne entière de 0
|
|
if (pivot_index < 0)
|
|
continue;
|
|
|
|
fullOrigin.at(pivot_index, 0) = origin.at(i, 0);
|
|
}
|
|
|
|
return {noyau, fullOrigin};
|
|
}
|
|
|
|
std::size_t Solver::Rank(Matrix&& a_Matrix) const {
|
|
return Image(std::move(a_Matrix)).GetCardinal();
|
|
}
|