fmatvec 0.0.0

linear_algebra_double.h
1 /* Copyright (C) 2003-2005 Martin Förg
2
3 * This library is free software; you can redistribute it and/or
4 * modify it under the terms of the GNU Lesser General Public
5 * License as published by the Free Software Foundation; either
6 * version 2.1 of the License, or (at your option) any later version.
7 *
8 * This library is distributed in the hope that it will be useful,
9 * but WITHOUT ANY WARRANTY; without even the implied warranty of
10 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11 * Lesser General Public License for more details.
12 *
13 * You should have received a copy of the GNU Lesser General Public
14 * License along with this library; if not, write to the Free Software
15 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
16 *
17 * Contact:
18 * martin.o.foerg@googlemail.com
19 *
20 */
21
22 #ifndef linear_algebrad_h
23 #define linear_algebrad_h
24
25 #include "square_matrix.h"
26 #include "vector.h"
27 #include "row_vector.h"
28 #include "diagonal_matrix.h"
29 #include "symmetric_matrix.h"
30 #include "band_matrix.h"
31 #include <complex>
32
33 //-------------------------------------
34 // Matrix operations
35 //-------------------------------------
36
37 namespace fmatvec {
38
42 FMATVEC_EXPORT int svd(Matrix<General, Ref, Ref, double> &A, Matrix<General, Ref, Ref, double> &S, SquareMatrix<Ref, double> &U, SquareMatrix<Ref, double> &VT, int Rueckgabe);
43
50 FMATVEC_EXPORT Vector<Ref, std::complex<double>> eigval(const SquareMatrix<Ref, double> &A);
51
60 FMATVEC_EXPORT int eigvec(const SquareMatrix<Ref, double> &A, SquareMatrix<Ref, std::complex<double>> &V, Vector<Ref, std::complex<double>> &w);
61
72 FMATVEC_EXPORT int eigvec(const Matrix<Symmetric, Ref, Ref, double> &A, const Matrix<Symmetric, Ref, Ref, double> &B, SquareMatrix<Ref, double> &eigenvectors, Vector<Ref, double> &eigenvalues);
73
80 FMATVEC_EXPORT Vector<Ref, double> eigval(const Matrix<Symmetric, Ref, Ref, double> &A);
81
93 FMATVEC_EXPORT Vector<Ref, double> eigvalSel(const Matrix<Symmetric, Ref, Ref, double> &A, int il, int iu, double abstol = 0);
94
116 FMATVEC_EXPORT int Doolittle_LU_with_Pivoting_Solve(const double *A, double B[], const int pivot[], double x[], int n);
117
127 FMATVEC_EXPORT Matrix<General, Ref, Ref, double> slvLU(const SquareMatrix<Ref, double> &A, const Matrix<General, Ref, Ref, double> &X);
128 FMATVEC_EXPORT Matrix<General, Var, Var, double> slvLU(const SquareMatrix<Var, double> &A, const Matrix<General, Var, Var, double> &X, int & info);
129
139 FMATVEC_EXPORT Vector<Ref, double> slvLU(const SquareMatrix<Ref, double> &A, const Vector<Ref, double> &x);
140
151 FMATVEC_EXPORT Vector<Ref, double> slvLU(const SquareMatrix<Ref, double> &A, const Vector<Ref, double> &x, int & info);
152 FMATVEC_EXPORT Vector<Var, double> slvLU(const SquareMatrix<Var, double> &A, const Vector<Var, double> &x, int & info);
153
164 FMATVEC_EXPORT Matrix<General, Ref, Ref, double> slvLUFac(const SquareMatrix<Ref, double> &A, const Matrix<General, Ref, Ref, double> &X, const Vector<Ref, int> &ipiv);
165 FMATVEC_EXPORT Matrix<General, Var, Var, double> slvLUFac(const SquareMatrix<Var, double> &A, const Matrix<General, Var, Var, double> &X, const Vector<Var, int> &ipiv);
166
177 FMATVEC_EXPORT Vector<Ref, double> slvLUFac(const SquareMatrix<Ref, double> &A, const Vector<Ref, double> &x, const Vector<Ref, int> &ipiv);
178 FMATVEC_EXPORT Vector<Var, double> slvLUFac(const SquareMatrix<Var, double> &A, const Vector<Var, double> &x, const Vector<Var, int> &ipiv);
179
189 FMATVEC_EXPORT Matrix<General, Ref, Ref, double> slvLL(const Matrix<Symmetric, Ref, Ref, double> &A, const Matrix<General, Ref, Ref, double> &X);
190
200 FMATVEC_EXPORT Vector<Ref, double> slvLL(const Matrix<Symmetric, Ref, Ref, double> &A, const Vector<Ref, double> &x);
201
211 FMATVEC_EXPORT Matrix<General, Ref, Ref, double> slvQR(const SquareMatrix<Ref, double> &A, const Matrix<General, Ref, Ref, double> &X);
212
222 FMATVEC_EXPORT Vector<Ref, double> slvQR(const SquareMatrix<Ref, double> &A, const Vector<Ref, double> &x);
223
233 FMATVEC_EXPORT Vector<Ref, double> slvQRLQ(const Matrix<General, Ref, Ref, double> &A, const Vector<Ref, double> &b);
234
244 FMATVEC_EXPORT Matrix<General, Ref, Ref, double> slvQRLQ(const Matrix<General, Ref, Ref, double> &A, const Matrix<General, Ref, Ref, double> &B);
245
253 FMATVEC_EXPORT SquareMatrix<Ref, double> inv(const SquareMatrix<Ref, double> &A);
254
262 FMATVEC_EXPORT Matrix<Symmetric, Ref, Ref, double> inv(const Matrix<Symmetric, Ref, Ref, double> &A);
263
271 FMATVEC_EXPORT Matrix<Diagonal, Ref, Ref, double> inv(const Matrix<Diagonal, Ref, Ref, double> &A);
272
281 FMATVEC_EXPORT Matrix<General, Ref, Ref, double> facLU(const Matrix<General, Ref, Ref, double> &A, Vector<Ref, int> &ipiv);
282
291 FMATVEC_EXPORT SquareMatrix<Ref, double> facLU(const SquareMatrix<Ref, double> &A, Vector<Ref, int> &ipiv);
292
293 FMATVEC_EXPORT SquareMatrix<Var, double> facLU(const SquareMatrix<Var, double> &A, Vector<Var, int> &ipiv);
294
306 FMATVEC_EXPORT int facLU(double *A, int pivot[], int n);
307
315 FMATVEC_EXPORT Matrix<Symmetric, Ref, Ref, double> facLL(const Matrix<Symmetric, Ref, Ref, double> &A);
316
323 FMATVEC_EXPORT double nrm1(const Vector<Ref, double> &x);
324
331 FMATVEC_EXPORT double nrm2(const Vector<Ref, double> &x);
332
339 FMATVEC_EXPORT double nrmInf(const Vector<Ref, double> &x);
340
348 FMATVEC_EXPORT double nrm1(const Matrix<General, Ref, Ref, double> &A);
349
356 FMATVEC_EXPORT double nrm2(const Matrix<General, Ref, Ref, double> &A);
357
365 FMATVEC_EXPORT double nrmInf(const Matrix<General, Ref, Ref, double> &A);
366
374 FMATVEC_EXPORT double nrmInf(const Matrix<Symmetric, Ref, Ref, double> &A);
375
382 FMATVEC_EXPORT double nrmFro(const Matrix<General, Ref, Ref, double> &A);
383
390 FMATVEC_EXPORT double rho(const SquareMatrix<Ref, double> &A);
391
398 FMATVEC_EXPORT double rho(const Matrix<Symmetric, Ref, Ref, double> &A);
399
400 FMATVEC_EXPORT Vector<Ref, double> slvLLFac(const Matrix<Symmetric, Ref, Ref, double> &A, const Vector<Ref, double> &x);
401
402 FMATVEC_EXPORT Matrix<General, Ref, Ref, double> slvLLFac(const Matrix<Symmetric, Ref, Ref, double> &A, const Matrix<General, Ref, Ref, double> &X);
403
404 FMATVEC_EXPORT Matrix<General, Ref, Ref, double> slvLS(const Matrix<General, Ref, Ref, double> &A, const Matrix<General, Ref, Ref, double> &B, double rcond = -1);
405
406 FMATVEC_EXPORT Vector<Ref, double> slvLS(const Matrix<General, Ref, Ref, double> &A, const Vector<Ref, double> &b, double rcond = -1);
407
408 //Matrix<General,Ref,Ref,double> slvLU(CBLAS_SIDE side, CBLAS_UPLO uplo, CBLAS_DIAG unit, const SquareMatrix<Ref,double> &A, const Matrix<General,Ref,Ref,double> &X, const Vector<Ref,int> &ipiv );
409
417 //Matrix<General,Ref,Ref,double> swap(const Matrix<General,Ref,Ref,double> &A, const Vector<Ref,int> &ipiv );
418
429 template <int size>
430 Vector<Fixed<size>, double> slvLU(const SquareMatrix<Fixed<size>, double> &A, const Vector<Fixed<size>, double> &b, int & info) {
431
432 if (size == 0)
433 return b;
434
435 int ipiv[size];
436
437 SquareMatrix<Fixed<size>, double> LU = A;
438
439 info = facLU(LU(), ipiv, size);
440
441 Vector<Fixed<size>, double> x;
442 Vector<Fixed<size>, double> y = b;
443
444 info = Doolittle_LU_with_Pivoting_Solve(LU(), y(), ipiv, x(), size);
445
446 return x;
447
448 }
449
450 template <int size>
451 SquareMatrix<Fixed<size>, double> facLU(const SquareMatrix<Fixed<size>, double> &A, Vector<Fixed<size>, int> &ipiv) {
452 SquareMatrix<Fixed<size>, double> LU = A;
453 facLU(LU(), ipiv(), size);
454 return LU;
455 }
456
457 template <int size>
458 Vector<Fixed<size>, double> slvLUFac(const SquareMatrix<Fixed<size>, double> &ALU, const Vector<Fixed<size>, double> &b, const Vector<Fixed<size>, int> &ipiv, int &info) {
459 Vector<Fixed<size>, double> x;
460 Vector<Fixed<size>, double> y = b;
461 info = Doolittle_LU_with_Pivoting_Solve(ALU(), y(), ipiv(), x(), size);
462 return x;
463 }
464
465 template <int size, int nrrhs>
466 Matrix<General, Fixed<size>, Fixed<nrrhs>, double> slvLUFac(const SquareMatrix<Fixed<size>, double> &ALU, const Matrix<General, Fixed<size>, Fixed<nrrhs>, double> &B, const Vector<Fixed<size>, int> &ipiv, int &info) {
467 Matrix<General, Fixed<size>, Fixed<nrrhs>, double> X;
468 Vector<Fixed<size>, double> x;
469 info = 0;
470 for(int c=0; c<nrrhs; c++) {
471 Vector<Fixed<size>, double> y = B.col(c);
472 int localInfo = Doolittle_LU_with_Pivoting_Solve(ALU(), y(), ipiv(), x(), size);
473 if(localInfo!=0)
474 info = localInfo;
475 X.set(c, x);
476 }
477 return X;
478 }
479
480 template <int size>
481 Matrix<General, Fixed<size>, Var, double> slvLUFac(const SquareMatrix<Fixed<size>, double> &ALU, const Matrix<General, Fixed<size>, Var, double> &B, const Vector<Fixed<size>, int> &ipiv, int &info) {
482 Matrix<General, Fixed<size>, Var, double> X(size,B.cols());
483 Vector<Fixed<size>, double> x;
484 info = 0;
485 for(int c=0; c<B.cols(); c++) {
486 Vector<Fixed<size>, double> y = B.col(c);
487 int localInfo = Doolittle_LU_with_Pivoting_Solve(ALU(), y(), ipiv(), x(), size);
488 if(localInfo!=0)
489 info = localInfo;
490 X.set(c, x);
491 }
492 return X;
493 }
494
495 }
496
497 #endif
fmatvec::svd
int svd(Matrix< General, Ref, Ref, double > &A, Matrix< General, Ref, Ref, double > &S, SquareMatrix< Ref, double > &U, SquareMatrix< Ref, double > &VT, int Rueckgabe)
Definition: linear_algebra_double.cc:779
fmatvec::facLL
Matrix< Symmetric, Ref, Ref, double > facLL(const Matrix< Symmetric, Ref, Ref, double > &A)
LL decomposition.
Definition: linear_algebra_double.cc:559
fmatvec::eigvec
int eigvec(const SquareMatrix< Ref, double > &A, SquareMatrix< Ref, std::complex< double >> &V, Vector< Ref, std::complex< double >> &w)
Eigenvectors and Eigenvalues.
Definition: linear_algebra_double.cc:715
fmatvec::Fixed
Definition: fmatvec.h:39
fmatvec::nrm1
double nrm1(const Vector< Ref, double > &x)
1-norm
Definition: linear_algebra_double.cc:578
fmatvec::eigvalSel
Vector< Ref, double > eigvalSel(const Matrix< Symmetric, Ref, Ref, double > &A, int il, int iu, double abstol)
Eigenvalues.
Definition: linear_algebra_double.cc:820
fmatvec::slvLU
Vector< Ref, std::complex< double > > slvLU(const SquareMatrix< Ref, std::complex< double >> &A, const Vector< Ref, std::complex< double >> &x)
System of linear equations.
Definition: linear_algebra_complex.cc:47
fmatvec::SquareMatrix
Definition: fmatvec.h:45
fmatvec::Vector
Definition: fmatvec.h:48
fmatvec::facLU
Matrix< General, Ref, Ref, double > facLU(const Matrix< General, Ref, Ref, double > &A, Vector< Ref, int > &ipiv)
LU decomposition.
Definition: linear_algebra_double.cc:441
fmatvec::Doolittle_LU_with_Pivoting_Solve
int Doolittle_LU_with_Pivoting_Solve(const double *A, double B[], const int pivot[], double x[], int n)
solve linear system with LU decomposed matrix
Definition: linear_algebra_double.cc:169
fmatvec::slvQR
Matrix< General, Ref, Ref, double > slvQR(const SquareMatrix< Ref, double > &A, const Matrix< General, Ref, Ref, double > &X)
Systems of linear equations.
Definition: linear_algebra_double.cc:150
fmatvec
Namespace fmatvec.
Definition: _memory.cc:28
fmatvec::slvLUFac
Matrix< General, Ref, Ref, double > slvLUFac(const SquareMatrix< Ref, double > &A, const Matrix< General, Ref, Ref, double > &X, const Vector< Ref, int > &ipiv)
Systems of linear equations.
Definition: linear_algebra_double.cc:103
fmatvec::eigval
Vector< Ref, std::complex< double > > eigval(const SquareMatrix< Ref, double > &A)
Eigenvalues.
Definition: linear_algebra_double.cc:695
fmatvec::inv
SquareMatrix< Ref, double > inv(const SquareMatrix< Ref, double > &A)
Inverse.
Definition: linear_algebra_double.cc:382
fmatvec::nrmFro
double nrmFro(const Matrix< General, Ref, Ref, double > &A)
Frobenius-norm.
Definition: linear_algebra_double.cc:892
fmatvec::slvLL
Vector< Ref, double > slvLL(const Matrix< Symmetric, Ref, Ref, double > &A, const Vector< Ref, double > &x)
System of linear equations.
Definition: linear_algebra_double.cc:603
fmatvec::slvQRLQ
Vector< Ref, double > slvQRLQ(const Matrix< General, Ref, Ref, double > &A, const Vector< Ref, double > &b)
System of overdetermined or underdetermined linear equations.
Definition: linear_algebra_double.cc:340
fmatvec::rho
double rho(const SquareMatrix< Ref, double > &A)
Spectral radius.
Definition: linear_algebra_double.cc:840