IT++ Logo
matfunc.cpp
Go to the documentation of this file.
1
30#include <itpp/base/matfunc.h>
33#include <limits>
34
36
37namespace itpp
38{
39
40// Square root of a square matrix (based on Octave sqrtm implementation)
41cmat sqrtm(const mat& A)
42{
43 return sqrtm(to_cmat(A));
44}
45
46// Square root of the complex square matrix A
47cmat sqrtm(const cmat& A)
48{
49 cmat U, T;
50 schur(A, U, T);
51
52 int n = U.rows();
53 cmat R(n, n);
54
55 R.zeros();
56 for (int j = 0; j < n; j++)
57 R(j, j) = std::sqrt(T(j, j));
58
59 const double fudge = std::sqrt(std::numeric_limits<double>::min());
60
61 for (int p = 0; p < n - 1; p++) {
62 for (int i = 0; i < n - (p + 1); i++) {
63 const int j = i + p + 1;
64 std::complex<double> s = T(i, j);
65
66 for (int k = i + 1; k < j; k++)
67 s -= R(i, k) * R(k, j);
68
69 const std::complex<double> d = R(i, i) + R(j, j) + fudge;
70 const std::complex<double> conj_d = conj(d);
71
72 R(i, j) = (s * conj_d) / (d * conj_d);
73 }
74 }
75
76 return U * R * U.H();
77}
78
79
80bool all(const bvec &testvec)
81{
82 for (int i = 0; i < testvec.length(); i++)
83 if (!testvec(i)) return false;
84 return true;
85}
86
87bool any(const bvec &testvec)
88{
89 for (int i = 0; i < testvec.length(); i++)
90 if (testvec(i)) return true;
91 return false;
92}
93
94// ----------------------------------------------------------------------
95// Instantiations
96// ----------------------------------------------------------------------
97
98template ITPP_EXPORT int length(const vec &v);
99template ITPP_EXPORT int length(const cvec &v);
100template ITPP_EXPORT int length(const svec &v);
101template ITPP_EXPORT int length(const ivec &v);
102template ITPP_EXPORT int length(const bvec &v);
103
104template ITPP_EXPORT double sum(const vec &v);
105template ITPP_EXPORT std::complex<double> sum(const cvec &v);
106template ITPP_EXPORT short sum(const svec &v);
107template ITPP_EXPORT int sum(const ivec &v);
108template ITPP_EXPORT bin sum(const bvec &v);
109
110template ITPP_EXPORT double sum_sqr(const vec &v);
111template ITPP_EXPORT std::complex<double> sum_sqr(const cvec &v);
112template ITPP_EXPORT short sum_sqr(const svec &v);
113template ITPP_EXPORT int sum_sqr(const ivec &v);
114template ITPP_EXPORT bin sum_sqr(const bvec &v);
115
116template ITPP_EXPORT vec cumsum(const vec &v);
117template ITPP_EXPORT cvec cumsum(const cvec &v);
118template ITPP_EXPORT svec cumsum(const svec &v);
119template ITPP_EXPORT ivec cumsum(const ivec &v);
120template ITPP_EXPORT bvec cumsum(const bvec &v);
121
122template ITPP_EXPORT double prod(const vec &v);
123template ITPP_EXPORT std::complex<double> prod(const cvec &v);
124template ITPP_EXPORT short prod(const svec &v);
125template ITPP_EXPORT int prod(const ivec &v);
126template ITPP_EXPORT bin prod(const bvec &v);
127
128template ITPP_EXPORT vec cross(const vec &v1, const vec &v2);
129template ITPP_EXPORT cvec cross(const cvec &v1, const cvec &v2);
130template ITPP_EXPORT ivec cross(const ivec &v1, const ivec &v2);
131template ITPP_EXPORT svec cross(const svec &v1, const svec &v2);
132template ITPP_EXPORT bvec cross(const bvec &v1, const bvec &v2);
133
134template ITPP_EXPORT vec reverse(const vec &in);
135template ITPP_EXPORT cvec reverse(const cvec &in);
136template ITPP_EXPORT svec reverse(const svec &in);
137template ITPP_EXPORT ivec reverse(const ivec &in);
138template ITPP_EXPORT bvec reverse(const bvec &in);
139
140template ITPP_EXPORT vec zero_pad(const vec &v, int n);
141template ITPP_EXPORT cvec zero_pad(const cvec &v, int n);
142template ITPP_EXPORT ivec zero_pad(const ivec &v, int n);
143template ITPP_EXPORT svec zero_pad(const svec &v, int n);
144template ITPP_EXPORT bvec zero_pad(const bvec &v, int n);
145
146template ITPP_EXPORT vec zero_pad(const vec &v);
147template ITPP_EXPORT cvec zero_pad(const cvec &v);
148template ITPP_EXPORT ivec zero_pad(const ivec &v);
149template ITPP_EXPORT svec zero_pad(const svec &v);
150template ITPP_EXPORT bvec zero_pad(const bvec &v);
151
152template ITPP_EXPORT mat zero_pad(const mat &, int, int);
153template ITPP_EXPORT cmat zero_pad(const cmat &, int, int);
154template ITPP_EXPORT imat zero_pad(const imat &, int, int);
155template ITPP_EXPORT smat zero_pad(const smat &, int, int);
156template ITPP_EXPORT bmat zero_pad(const bmat &, int, int);
157
158template ITPP_EXPORT vec sum(const mat &m, int dim);
159template ITPP_EXPORT cvec sum(const cmat &m, int dim);
160template ITPP_EXPORT svec sum(const smat &m, int dim);
161template ITPP_EXPORT ivec sum(const imat &m, int dim);
162template ITPP_EXPORT bvec sum(const bmat &m, int dim);
163
164template ITPP_EXPORT double sumsum(const mat &X);
165template ITPP_EXPORT std::complex<double> sumsum(const cmat &X);
166template ITPP_EXPORT short sumsum(const smat &X);
167template ITPP_EXPORT int sumsum(const imat &X);
168template ITPP_EXPORT bin sumsum(const bmat &X);
169
170template ITPP_EXPORT vec sum_sqr(const mat & m, int dim);
171template ITPP_EXPORT cvec sum_sqr(const cmat &m, int dim);
172template ITPP_EXPORT svec sum_sqr(const smat &m, int dim);
173template ITPP_EXPORT ivec sum_sqr(const imat &m, int dim);
174template ITPP_EXPORT bvec sum_sqr(const bmat &m, int dim);
175
176template ITPP_EXPORT mat cumsum(const mat &m, int dim);
177template ITPP_EXPORT cmat cumsum(const cmat &m, int dim);
178template ITPP_EXPORT smat cumsum(const smat &m, int dim);
179template ITPP_EXPORT imat cumsum(const imat &m, int dim);
180template ITPP_EXPORT bmat cumsum(const bmat &m, int dim);
181
182template ITPP_EXPORT vec prod(const mat &m, int dim);
183template ITPP_EXPORT cvec prod(const cmat &v, int dim);
184template ITPP_EXPORT svec prod(const smat &m, int dim);
185template ITPP_EXPORT ivec prod(const imat &m, int dim);
186template ITPP_EXPORT bvec prod(const bmat &m, int dim);
187
188template ITPP_EXPORT vec diag(const mat &in);
189template ITPP_EXPORT cvec diag(const cmat &in);
190
191template ITPP_EXPORT void diag(const vec &in, mat &m);
192template ITPP_EXPORT void diag(const cvec &in, cmat &m);
193
194template ITPP_EXPORT mat diag(const vec &v, const int K);
195template ITPP_EXPORT cmat diag(const cvec &v, const int K);
196
197template ITPP_EXPORT mat bidiag(const vec &, const vec &);
198template ITPP_EXPORT cmat bidiag(const cvec &, const cvec &);
199
200template ITPP_EXPORT void bidiag(const vec &, const vec &, mat &);
201template ITPP_EXPORT void bidiag(const cvec &, const cvec &, cmat &);
202
203template ITPP_EXPORT void bidiag(const mat &, vec &, vec &);
204template ITPP_EXPORT void bidiag(const cmat &, cvec &, cvec &);
205
206template ITPP_EXPORT mat tridiag(const vec &main, const vec &, const vec &);
207template ITPP_EXPORT cmat tridiag(const cvec &main, const cvec &, const cvec &);
208
209template ITPP_EXPORT void tridiag(const vec &main, const vec &, const vec &, mat &);
210template ITPP_EXPORT void tridiag(const cvec &main, const cvec &, const cvec &, cmat &);
211
212template ITPP_EXPORT void tridiag(const mat &m, vec &, vec &, vec &);
213template ITPP_EXPORT void tridiag(const cmat &m, cvec &, cvec &, cvec &);
214
215template ITPP_EXPORT double trace(const mat &in);
216template ITPP_EXPORT std::complex<double> trace(const cmat &in);
217template ITPP_EXPORT short trace(const smat &in);
218template ITPP_EXPORT int trace(const imat &in);
219template ITPP_EXPORT bin trace(const bmat &in);
220
221template ITPP_EXPORT void transpose(const mat &m, mat &out);
222template ITPP_EXPORT void transpose(const cmat &m, cmat &out);
223template ITPP_EXPORT void transpose(const smat &m, smat &out);
224template ITPP_EXPORT void transpose(const imat &m, imat &out);
225template ITPP_EXPORT void transpose(const bmat &m, bmat &out);
226
227template ITPP_EXPORT mat transpose(const mat &m);
228template ITPP_EXPORT cmat transpose(const cmat &m);
229template ITPP_EXPORT smat transpose(const smat &m);
230template ITPP_EXPORT imat transpose(const imat &m);
231template ITPP_EXPORT bmat transpose(const bmat &m);
232
233template ITPP_EXPORT void hermitian_transpose(const mat &m, mat &out);
234template ITPP_EXPORT void hermitian_transpose(const cmat &m, cmat &out);
235template ITPP_EXPORT void hermitian_transpose(const smat &m, smat &out);
236template ITPP_EXPORT void hermitian_transpose(const imat &m, imat &out);
237template ITPP_EXPORT void hermitian_transpose(const bmat &m, bmat &out);
238
239template ITPP_EXPORT mat hermitian_transpose(const mat &m);
240template ITPP_EXPORT cmat hermitian_transpose(const cmat &m);
241template ITPP_EXPORT smat hermitian_transpose(const smat &m);
242template ITPP_EXPORT imat hermitian_transpose(const imat &m);
243template ITPP_EXPORT bmat hermitian_transpose(const bmat &m);
244
245template ITPP_EXPORT bool is_hermitian(const mat &X);
246template ITPP_EXPORT bool is_hermitian(const cmat &X);
247
248template ITPP_EXPORT bool is_unitary(const mat &X);
249template ITPP_EXPORT bool is_unitary(const cmat &X);
250
251template ITPP_EXPORT vec rvectorize(const mat &m);
252template ITPP_EXPORT cvec rvectorize(const cmat &m);
253template ITPP_EXPORT ivec rvectorize(const imat &m);
254template ITPP_EXPORT svec rvectorize(const smat &m);
255template ITPP_EXPORT bvec rvectorize(const bmat &m);
256
257template ITPP_EXPORT vec cvectorize(const mat &m);
258template ITPP_EXPORT cvec cvectorize(const cmat &m);
259template ITPP_EXPORT ivec cvectorize(const imat &m);
260template ITPP_EXPORT svec cvectorize(const smat &m);
261template ITPP_EXPORT bvec cvectorize(const bmat &m);
262
263template ITPP_EXPORT mat reshape(const mat &m, int rows, int cols);
264template ITPP_EXPORT cmat reshape(const cmat &m, int rows, int cols);
265template ITPP_EXPORT imat reshape(const imat &m, int rows, int cols);
266template ITPP_EXPORT smat reshape(const smat &m, int rows, int cols);
267template ITPP_EXPORT bmat reshape(const bmat &m, int rows, int cols);
268
269template ITPP_EXPORT mat reshape(const vec &m, int rows, int cols);
270template ITPP_EXPORT cmat reshape(const cvec &m, int rows, int cols);
271template ITPP_EXPORT imat reshape(const ivec &m, int rows, int cols);
272template ITPP_EXPORT smat reshape(const svec &m, int rows, int cols);
273template ITPP_EXPORT bmat reshape(const bvec &m, int rows, int cols);
274
275template ITPP_EXPORT mat kron(const mat &X, const mat &Y);
276template ITPP_EXPORT cmat kron(const cmat &X, const cmat &Y);
277template ITPP_EXPORT imat kron(const imat &X, const imat &Y);
278template ITPP_EXPORT smat kron(const smat &X, const smat &Y);
279template ITPP_EXPORT bmat kron(const bmat &X, const bmat &Y);
280
281template ITPP_EXPORT vec repmat(const vec &v, int n);
282template ITPP_EXPORT cvec repmat(const cvec &v, int n);
283template ITPP_EXPORT ivec repmat(const ivec &v, int n);
284template ITPP_EXPORT svec repmat(const svec &v, int n);
285template ITPP_EXPORT bvec repmat(const bvec &v, int n);
286
287template ITPP_EXPORT mat repmat(const vec &v, int m, int n, bool transpose);
288template ITPP_EXPORT cmat repmat(const cvec &v, int m, int n, bool transpose);
289template ITPP_EXPORT imat repmat(const ivec &v, int m, int n, bool transpose);
290template ITPP_EXPORT smat repmat(const svec &v, int m, int n, bool transpose);
291template ITPP_EXPORT bmat repmat(const bvec &v, int m, int n, bool transpose);
292
293template ITPP_EXPORT mat repmat(const mat &data, int m, int n);
294template ITPP_EXPORT cmat repmat(const cmat &data, int m, int n);
295template ITPP_EXPORT imat repmat(const imat &data, int m, int n);
296template ITPP_EXPORT smat repmat(const smat &data, int m, int n);
297template ITPP_EXPORT bmat repmat(const bmat &data, int m, int n);
298
299} // namespace itpp
300
Definitions of converters between different vector and matrix types.
T trace(const Mat< T > &m)
The trace of the matrix m, i.e. the sum of the diagonal elements.
Definition matfunc.h:762
Mat< T > tridiag(const Vec< T > &main, const Vec< T > &sup, const Vec< T > &sub)
Returns a matrix with the elements of main on the diagonal, the elements of sup on the diagonal row a...
Definition matfunc.h:690
Mat< T > diag(const Vec< T > &v, const int K=0)
Create a diagonal matrix using vector v as its diagonal.
Definition matfunc.h:557
Mat< T > bidiag(const Vec< T > &main, const Vec< T > &sup)
Returns a matrix with the elements of the input vector main on the diagonal and the elements of the i...
Definition matfunc.h:617
Mat< Num_T > kron(const Mat< Num_T > &X, const Mat< Num_T > &Y)
Computes the Kronecker product of two matrices.
Definition matfunc.h:443
T sumsum(const Mat< T > &X)
Sum of all elements in the given matrix. Fast version of sum(sum(X))
Definition matfunc.h:101
ITPP_EXPORT cmat sqrtm(const cmat &A)
Square root of the complex square matrix A.
bool is_unitary(const Mat< Num_T > &X)
Returns true if matrix X is unitary, false otherwise.
Definition matfunc.h:355
Vec< T > zero_pad(const Vec< T > &v, int n)
Zero-pad a vector to size n.
Definition matfunc.h:257
Vec< T > cumsum(const Vec< T > &v)
Cumulative sum of all elements in the vector.
Definition matfunc.h:157
T sum(const Vec< T > &v)
Sum of all elements in the vector.
Definition matfunc.h:59
bool is_hermitian(const Mat< Num_T > &X)
Returns true if matrix X is hermitian, false otherwise.
Definition matfunc.h:336
void transpose(const Mat< T > &m, Mat< T > &out)
Transposition of the matrix m returning the transposed matrix in out.
Definition matfunc.h:308
T sum_sqr(const Vec< T > &v)
Sum of square of the elements in a vector.
Definition matfunc.h:116
int length(const Vec< T > &v)
Length of vector.
Definition matfunc.h:51
void hermitian_transpose(const Mat< T > &m, Mat< T > &out)
Definition matfunc.h:318
T prod(const Vec< T > &v)
The product of all elements in the vector.
Definition matfunc.h:195
Vec< T > cross(const Vec< T > &v1, const Vec< T > &v2)
Vector cross product. Vectors need to be of size 3.
Definition matfunc.h:240
bool schur(const mat &A, mat &U, mat &T)
Schur decomposition of a real matrix.
Definition schur.cpp:104
cvec conj(const cvec &x)
Conjugate of complex value.
Definition elem_math.h:226
Vec< T > rvectorize(const Mat< T > &m)
Row vectorize the matrix [(0,0) (0,1) ... (N-1,N-2) (N-1,N-1)].
Definition matfunc.h:789
Vec< T > cvectorize(const Mat< T > &m)
Column vectorize the matrix [(0,0) (1,0) ... (N-2,N-1) (N-1,N-1)].
Definition matfunc.h:803
Vec< T > reverse(const Vec< T > &in)
Reverse the input vector.
Definition matfunc.h:777
Mat< T > reshape(const Mat< T > &m, int rows, int cols)
Reshape the matrix into an rows*cols matrix.
Definition matfunc.h:822
Mat< bin > bmat
bin matrix
Definition mat.h:508
Various functions on vectors and matrices - header file.
itpp namespace
Definition itmex.h:37
cmat to_cmat(const Mat< T > &m)
Converts a Mat<T> to cmat.
Definition converters.h:232
ITPP_EXPORT bool any(const bvec &testvec)
Returns true if any element is one and false otherwise.
ITPP_EXPORT bool all(const bvec &testvec)
Returns true if all elements are ones and false otherwise.
Definitions of Schur decomposition functions.

Generated on Tue Aug 17 2021 10:59:15 for IT++ by Doxygen 1.9.8