30 # include <itpp/config.h> 32 # include <itpp/config_msvc.h> 35 #if defined(HAVE_LAPACK) 36 # include <itpp/base/algebra/lapack.h> 46 #if defined(HAVE_LAPACK) 48 bool eig_sym(
const mat &A, vec &d, mat &V)
50 it_assert_debug(A.rows() == A.cols(),
"eig_sym: Matrix is not symmetric");
54 char jobz =
'V', uplo =
'U';
55 int n, lda, lwork, info;
64 dsyev_(&jobz, &uplo, &n, V._data(), &lda, d._data(), work._data(), &lwork, &info);
69 bool eig_sym(
const mat &A, vec &d)
71 it_assert_debug(A.rows() == A.cols(),
"eig_sym: Matrix is not symmetric");
75 char jobz =
'N', uplo =
'U';
76 int n, lda, lwork, info;
85 dsyev_(&jobz, &uplo, &n, B._data(), &lda, d._data(), work._data(), &lwork, &info);
90 bool eig_sym(
const cmat &A, vec &d, cmat &V)
92 it_assert_debug(A.rows() == A.cols(),
"eig_sym: Matrix is not hermitian");
96 char jobz =
'V', uplo =
'U';
97 int n, lda, lwork, info;
101 d.set_size(n,
false);
107 zheev_(&jobz, &uplo, &n, V._data(), &lda, d._data(), work._data(), &lwork, rwork._data(), &info);
112 bool eig_sym(
const cmat &A, vec &d)
114 it_assert_debug(A.rows() == A.cols(),
"eig_sym: Matrix is not hermitian");
118 char jobz =
'N', uplo =
'U';
119 int n, lda, lwork, info;
123 d.set_size(n,
false);
129 zheev_(&jobz, &uplo, &n, B._data(), &lda, d._data(), work._data(), &lwork, rwork._data(), &info);
136 bool eig(
const mat &A, cvec &d, cmat &V)
140 char jobvl =
'N', jobvr =
'V';
141 int n, lda, ldvl, ldvr, lwork, info;
154 dgeev_(&jobvl, &jobvr, &n, B._data(), &lda, wr._data(), wi._data(), vl._data(), &ldvl, vr._data(), &ldvr, work._data(), &lwork, &info);
159 V.set_size(n, n,
false);
160 for (
int j = 0; j < n; j++) {
162 if ((j < n - 1) && d(j) ==
std::conj(d(j + 1))) {
163 V.set_col(j,
to_cvec(vr.get_col(j), vr.get_col(j + 1)));
164 V.set_col(j + 1,
to_cvec(vr.get_col(j), -vr.get_col(j + 1)));
168 V.set_col(j,
to_cvec(vr.get_col(j)));
176 bool eig(
const mat &A, cvec &d)
180 char jobvl =
'N', jobvr =
'N';
181 int n, lda, ldvl, ldvr, lwork, info;
194 dgeev_(&jobvl, &jobvr, &n, B._data(), &lda, wr._data(), wi._data(), vl._data(), &ldvl, vr._data(), &ldvr, work._data(), &lwork, &info);
201 bool eig(
const cmat &A, cvec &d, cmat &V)
205 char jobvl =
'N', jobvr =
'V';
206 int n, lda, ldvl, ldvr, lwork, info;
212 d.set_size(n,
false);
213 V.set_size(n, n,
false);
220 zgeev_(&jobvl, &jobvr, &n, B._data(), &lda, d._data(), vl._data(), &ldvl, V._data(), &ldvr, work._data(), &lwork, rwork._data(), &info);
226 bool eig(
const cmat &A, cvec &d)
230 char jobvl =
'N', jobvr =
'N';
231 int n, lda, ldvl, ldvr, lwork, info;
237 d.set_size(n,
false);
244 zgeev_(&jobvl, &jobvr, &n, B._data(), &lda, d._data(), vl._data(), &ldvl, vr._data(), &ldvr, work._data(), &lwork, rwork._data(), &info);
254 it_error(
"LAPACK library is needed to use eig_sym() function");
260 it_error(
"LAPACK library is needed to use eig_sym() function");
266 it_error(
"LAPACK library is needed to use eig_sym() function");
272 it_error(
"LAPACK library is needed to use eig_sym() function");
277 bool eig(
const mat &A, cvec &d, cmat &V)
279 it_error(
"LAPACK library is needed to use eig() function");
283 bool eig(
const mat &A, cvec &d)
285 it_error(
"LAPACK library is needed to use eig() function");
289 bool eig(
const cmat &A, cvec &d, cmat &V)
291 it_error(
"LAPACK library is needed to use eig() function");
295 bool eig(
const cmat &A, cvec &d)
297 it_error(
"LAPACK library is needed to use eig() function");
301 #endif // HAVE_LAPACK Definitions of eigenvalue decomposition functions.
cvec conj(const cvec &x)
Conjugate of complex value.
#define it_assert_debug(t, s)
Abort if t is not true and NDEBUG is not defined.
Definitions of converters between different vector and matrix types.
T max(const Vec< T > &v)
Maximum value of vector.
bool eig(const mat &A, cvec &d, cmat &V)
Calculates the eigenvalues and eigenvectors of a real non-symmetric matrix.
bool eig_sym(const mat &A, vec &d, mat &V)
Calculates the eigenvalues and eigenvectors of a symmetric real matrix.
#define it_error(s)
Abort unconditionally.
cvec to_cvec(const Vec< T > &v)
Converts a Vec<T> to cvec.