Go to the documentation of this file.
26 #ifndef O2SCL_LINEAR_SOLVER_H
27 #define O2SCL_LINEAR_SOLVER_H
29 #include <gsl/gsl_linalg.h>
31 #include <o2scl/permutation.h>
51 template<
class vec_t=boost::numeric::ublas::vector<
double>,
52 class mat_t=boost::numeric::ublas::matrix<
double> >
60 virtual void solve(
size_t n, mat_t &a, vec_t &b, vec_t &x)=0;
65 template <
class vec_t=boost::numeric::ublas::vector<
double>,
66 class mat_t=boost::numeric::ublas::matrix<
double> >
72 virtual void solve(
size_t n, mat_t &A,
87 template <
class vec_t=boost::numeric::ublas::vector<
double>,
88 class mat_t=boost::numeric::ublas::matrix<
double> >
94 virtual void solve(
size_t n, mat_t &A, vec_t &b, vec_t &x) {
107 template <
class vec_t=boost::numeric::ublas::vector<
double>,
108 class mat_t=boost::numeric::ublas::matrix<
double> >
115 virtual void solve(
size_t n, mat_t &A, vec_t &b, vec_t &x) {
127 #if defined (O2SCL_COND_FLAG) || defined (DOXYGEN)
128 #if defined (O2SCL_ARMA) || defined (DOXYGEN)
138 virtual void solve(
size_t n, arma_mat_t &A, arma_vec_t &b,
147 #if defined (O2SCL_EIGEN) || defined (DOXYGEN)
148 #include <eigen3/Eigen/Dense>
157 template<
class eigen_vec_t,
class eigen_mat_t>
160 virtual void solve(
size_t n, eigen_mat_t &A, eigen_vec_t &b,
162 x=A.householderQr().solve(b);
174 template<
class eigen_vec_t,
class eigen_mat_t>
177 virtual void solve(
size_t n, eigen_mat_t &A, eigen_vec_t &b,
179 x=A.colPivHouseholderQr().solve(b);
191 template<
class eigen_vec_t,
class eigen_mat_t>
194 virtual void solve(
size_t n, eigen_mat_t &A, eigen_vec_t &b,
196 x=A.fullPivHouseholderQr().solve(b);
210 template<
class eigen_vec_t,
class eigen_mat_t>
213 virtual void solve(
size_t n, eigen_mat_t &A, eigen_vec_t &b,
215 x=A.partialPivLu().solve(b);
227 template<
class eigen_vec_t,
class eigen_mat_t>
230 virtual void solve(
size_t n, eigen_mat_t &A, eigen_vec_t &b,
232 x=A.fullPivLu().solve(b);
246 template<
class eigen_vec_t,
class eigen_mat_t>
249 virtual void solve(
size_t n, eigen_mat_t &A, eigen_vec_t &b,
265 template<
class eigen_vec_t,
class eigen_mat_t>
268 virtual void solve(
size_t n, eigen_mat_t &A, eigen_vec_t &b,
278 #include <o2scl/linear_special.h>
virtual void solve(size_t n, mat_t &A, vec_t &b, vec_t &x)
Solve square linear system of size n.
virtual void solve(size_t n, eigen_mat_t &A, eigen_vec_t &b, eigen_vec_t &x)
Solve square linear system of size n.
virtual void solve(size_t n, eigen_mat_t &A, eigen_vec_t &b, eigen_vec_t &x)
Solve square linear system of size n.
virtual void solve(size_t n, mat_t &A, vec_t &b, vec_t &x)
Solve square linear system of size n.
The namespace for linear algebra classes and functions.
A class for representing permutations.
Eigen linear solver using QR decomposition with column pivoting.
A generic solver for the linear system [abstract base].
Generic linear solver using QR decomposition.
virtual void solve(size_t n, eigen_mat_t &A, eigen_vec_t &b, eigen_vec_t &x)
Solve square linear system of size n.
int LU_decomp(const size_t N, mat_t &A, o2scl::permutation &p, int &signum)
Compute the LU decomposition of the matrix A.
virtual void solve(size_t n, eigen_mat_t &A, eigen_vec_t &b, eigen_vec_t &x)
Solve square linear system of size n.
Eigen linear solver using LLT decomposition with full pivoting.
virtual void solve(size_t n, eigen_mat_t &A, eigen_vec_t &b, eigen_vec_t &x)
Solve square linear system of size n.
void QR_decomp(size_t M, size_t N, mat_t &A, vec_t &tau)
Compute the QR decomposition of matrix A.
virtual void solve(size_t n, mat_t &a, vec_t &b, vec_t &x)=0
Solve square linear system of size n.
Eigen linear solver using LDLT decomposition with full pivoting.
Generic Householder linear solver.
virtual void solve(size_t n, eigen_mat_t &A, eigen_vec_t &b, eigen_vec_t &x)
Solve square linear system of size n.
int HH_solve(size_t n, mat_t &A, const vec_t &b, vec2_t &x)
Solve linear system after Householder decomposition.
virtual void solve(size_t n, eigen_mat_t &A, eigen_vec_t &b, eigen_vec_t &x)
Solve square linear system of size n.
Eigen linear solver using QR decomposition with full pivoting.
Eigen linear solver using QR decomposition with column pivoting.
Generic linear solver using LU decomposition.
virtual void solve(size_t n, arma_mat_t &A, arma_vec_t &b, arma_vec_t &x)
Solve square linear system of size n.
int LU_solve(const size_t N, const mat_t &LU, const o2scl::permutation &p, const vec_t &b, vec2_t &x)
Solve a linear system after LU decomposition.
virtual void solve(size_t n, mat_t &A, vec_t &b, vec_t &x)
Solve square linear system of size n.
void QR_solve(size_t N, const mat_t &QR, const vec_t &tau, const vec2_t &b, vec3_t &x)
Solve the system A x = b using the QR factorization.
Eigen linear solver using LU decomposition with full pivoting.
Eigen linear solver using LU decomposition with partial pivoting.
Documentation generated with Doxygen. Provided under the
GNU Free Documentation License (see License Information).