Linear Algebra with C++

[VeganFanatic]

I have been working on a project for the last year and I have make much progress but much remains. First I will provide headers of what I have done so far.

template <class T> class vector; // vector space
template <class T> class matrix; // derived from vector
template <typename T> T norm(vector<T> v);
template <typename T> T norm(matrix<T> m);
template <typename T> void transpose(matrix<T> &m);
template <typename T> void conjugate_transpose(matrix<T> &m);
template <typename T> bool hermetian(matrix<T> &a, matrix<T> &b);
template <typename T> matrix<T> tensor_product(vector<T> row, vector<T> col);
template <typename T> matrix<T> minor(matrix<T> &source, int row, int col);
template <typename T> bool diagonal_dominant(matrix<T> a);
template <typename T> matrix<T> gauss_elimination(matrix<T> &m);
template <typename T> matrix<T> gauss_jordan(matrix<T> &m);
template <typename T> void gauss_siedel(matrix<T> a, vector<T> b, vector<T> x, double tolerance);

Now I can use vector and use it with math equations fine. Same for the matrix template. So given my work to make these helpful templates, I am now out of my math book big time and into the realm of advanced linear algebra. This is my Gauss elimination, and I was looking to leverage it to do LU decomposition etc to solve Ax = b more expediently.

template <typename T> matrix<T> gauss_elimination(matrix<T> &mat) {
int i = 0; int j = 0; int maxi;
matrix<T> a = mat;
int m = (int)a.size();
int n = (int)a[0].size();
while ((i < m) && (j < n)) {
maxi = i;
// sweep down the column to find the pivot element
for (int k = i + 1; k < m; k++) // partial pivoting
if (std::abs(a[k][j]) > std::abs(a[maxi][j])) maxi = k;
// eliminate
if (a[maxi][j] != 0) {
swap(a[i], a[maxi]); // row op 1, swap
a[i] = a[i] / a[i][j]; // row op 2, scale the pivot to 1
for (int u = i + 1; u < m; u++)
a[u] = a[u] - (a[i] * a[u][j]); // row op 3, zero out the column below
i++;
}
j++;
}
return a;
}

Suggestions?

http://www.contract-developer.tk

[jean Davy]

My suggestion is to not reinvent the wheel... ;) You should have a look at boost Basic Linear Algebra Library[^] which is a part of Boost.Math[^] I like also Generic Geometry Library[^] Regards,