IT++ Logo
matfunc.cpp
Go to the documentation of this file.
1 
30 #include <itpp/base/matfunc.h>
31 #include <itpp/base/algebra/schur.h>
32 #include <itpp/base/converters.h>
33 #include <limits>
34 
36 
37 namespace itpp
38 {
39 
40 // Square root of a square matrix (based on Octave sqrtm implementation)
41 cmat sqrtm(const mat& A)
42 {
43  return sqrtm(to_cmat(A));
44 }
45 
46 // Square root of the complex square matrix A
47 cmat 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 
80 bool 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 
87 bool 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 
98 template ITPP_EXPORT int length(const vec &v);
99 template ITPP_EXPORT int length(const cvec &v);
100 template ITPP_EXPORT int length(const svec &v);
101 template ITPP_EXPORT int length(const ivec &v);
102 template ITPP_EXPORT int length(const bvec &v);
103 
104 template ITPP_EXPORT double sum(const vec &v);
105 template ITPP_EXPORT std::complex<double> sum(const cvec &v);
106 template ITPP_EXPORT short sum(const svec &v);
107 template ITPP_EXPORT int sum(const ivec &v);
108 template ITPP_EXPORT bin sum(const bvec &v);
109 
110 template ITPP_EXPORT double sum_sqr(const vec &v);
111 template ITPP_EXPORT std::complex<double> sum_sqr(const cvec &v);
112 template ITPP_EXPORT short sum_sqr(const svec &v);
113 template ITPP_EXPORT int sum_sqr(const ivec &v);
114 template ITPP_EXPORT bin sum_sqr(const bvec &v);
115 
116 template ITPP_EXPORT vec cumsum(const vec &v);
117 template ITPP_EXPORT cvec cumsum(const cvec &v);
118 template ITPP_EXPORT svec cumsum(const svec &v);
119 template ITPP_EXPORT ivec cumsum(const ivec &v);
120 template ITPP_EXPORT bvec cumsum(const bvec &v);
121 
122 template ITPP_EXPORT double prod(const vec &v);
123 template ITPP_EXPORT std::complex<double> prod(const cvec &v);
124 template ITPP_EXPORT short prod(const svec &v);
125 template ITPP_EXPORT int prod(const ivec &v);
126 template ITPP_EXPORT bin prod(const bvec &v);
127 
128 template ITPP_EXPORT vec cross(const vec &v1, const vec &v2);
129 template ITPP_EXPORT cvec cross(const cvec &v1, const cvec &v2);
130 template ITPP_EXPORT ivec cross(const ivec &v1, const ivec &v2);
131 template ITPP_EXPORT svec cross(const svec &v1, const svec &v2);
132 template ITPP_EXPORT bvec cross(const bvec &v1, const bvec &v2);
133 
134 template ITPP_EXPORT vec reverse(const vec &in);
135 template ITPP_EXPORT cvec reverse(const cvec &in);
136 template ITPP_EXPORT svec reverse(const svec &in);
137 template ITPP_EXPORT ivec reverse(const ivec &in);
138 template ITPP_EXPORT bvec reverse(const bvec &in);
139 
140 template ITPP_EXPORT vec zero_pad(const vec &v, int n);
141 template ITPP_EXPORT cvec zero_pad(const cvec &v, int n);
142 template ITPP_EXPORT ivec zero_pad(const ivec &v, int n);
143 template ITPP_EXPORT svec zero_pad(const svec &v, int n);
144 template ITPP_EXPORT bvec zero_pad(const bvec &v, int n);
145 
146 template ITPP_EXPORT vec zero_pad(const vec &v);
147 template ITPP_EXPORT cvec zero_pad(const cvec &v);
148 template ITPP_EXPORT ivec zero_pad(const ivec &v);
149 template ITPP_EXPORT svec zero_pad(const svec &v);
150 template ITPP_EXPORT bvec zero_pad(const bvec &v);
151 
152 template ITPP_EXPORT mat zero_pad(const mat &, int, int);
153 template ITPP_EXPORT cmat zero_pad(const cmat &, int, int);
154 template ITPP_EXPORT imat zero_pad(const imat &, int, int);
155 template ITPP_EXPORT smat zero_pad(const smat &, int, int);
156 template ITPP_EXPORT bmat zero_pad(const bmat &, int, int);
157 
158 template ITPP_EXPORT vec sum(const mat &m, int dim);
159 template ITPP_EXPORT cvec sum(const cmat &m, int dim);
160 template ITPP_EXPORT svec sum(const smat &m, int dim);
161 template ITPP_EXPORT ivec sum(const imat &m, int dim);
162 template ITPP_EXPORT bvec sum(const bmat &m, int dim);
163 
164 template ITPP_EXPORT double sumsum(const mat &X);
165 template ITPP_EXPORT std::complex<double> sumsum(const cmat &X);
166 template ITPP_EXPORT short sumsum(const smat &X);
167 template ITPP_EXPORT int sumsum(const imat &X);
168 template ITPP_EXPORT bin sumsum(const bmat &X);
169 
170 template ITPP_EXPORT vec sum_sqr(const mat & m, int dim);
171 template ITPP_EXPORT cvec sum_sqr(const cmat &m, int dim);
172 template ITPP_EXPORT svec sum_sqr(const smat &m, int dim);
173 template ITPP_EXPORT ivec sum_sqr(const imat &m, int dim);
174 template ITPP_EXPORT bvec sum_sqr(const bmat &m, int dim);
175 
176 template ITPP_EXPORT mat cumsum(const mat &m, int dim);
177 template ITPP_EXPORT cmat cumsum(const cmat &m, int dim);
178 template ITPP_EXPORT smat cumsum(const smat &m, int dim);
179 template ITPP_EXPORT imat cumsum(const imat &m, int dim);
180 template ITPP_EXPORT bmat cumsum(const bmat &m, int dim);
181 
182 template ITPP_EXPORT vec prod(const mat &m, int dim);
183 template ITPP_EXPORT cvec prod(const cmat &v, int dim);
184 template ITPP_EXPORT svec prod(const smat &m, int dim);
185 template ITPP_EXPORT ivec prod(const imat &m, int dim);
186 template ITPP_EXPORT bvec prod(const bmat &m, int dim);
187 
188 template ITPP_EXPORT vec diag(const mat &in);
189 template ITPP_EXPORT cvec diag(const cmat &in);
190 
191 template ITPP_EXPORT void diag(const vec &in, mat &m);
192 template ITPP_EXPORT void diag(const cvec &in, cmat &m);
193 
194 template ITPP_EXPORT mat diag(const vec &v, const int K);
195 template ITPP_EXPORT cmat diag(const cvec &v, const int K);
196 
197 template ITPP_EXPORT mat bidiag(const vec &, const vec &);
198 template ITPP_EXPORT cmat bidiag(const cvec &, const cvec &);
199 
200 template ITPP_EXPORT void bidiag(const vec &, const vec &, mat &);
201 template ITPP_EXPORT void bidiag(const cvec &, const cvec &, cmat &);
202 
203 template ITPP_EXPORT void bidiag(const mat &, vec &, vec &);
204 template ITPP_EXPORT void bidiag(const cmat &, cvec &, cvec &);
205 
206 template ITPP_EXPORT mat tridiag(const vec &main, const vec &, const vec &);
207 template ITPP_EXPORT cmat tridiag(const cvec &main, const cvec &, const cvec &);
208 
209 template ITPP_EXPORT void tridiag(const vec &main, const vec &, const vec &, mat &);
210 template ITPP_EXPORT void tridiag(const cvec &main, const cvec &, const cvec &, cmat &);
211 
212 template ITPP_EXPORT void tridiag(const mat &m, vec &, vec &, vec &);
213 template ITPP_EXPORT void tridiag(const cmat &m, cvec &, cvec &, cvec &);
214 
215 template ITPP_EXPORT double trace(const mat &in);
216 template ITPP_EXPORT std::complex<double> trace(const cmat &in);
217 template ITPP_EXPORT short trace(const smat &in);
218 template ITPP_EXPORT int trace(const imat &in);
219 template ITPP_EXPORT bin trace(const bmat &in);
220 
221 template ITPP_EXPORT void transpose(const mat &m, mat &out);
222 template ITPP_EXPORT void transpose(const cmat &m, cmat &out);
223 template ITPP_EXPORT void transpose(const smat &m, smat &out);
224 template ITPP_EXPORT void transpose(const imat &m, imat &out);
225 template ITPP_EXPORT void transpose(const bmat &m, bmat &out);
226 
227 template ITPP_EXPORT mat transpose(const mat &m);
228 template ITPP_EXPORT cmat transpose(const cmat &m);
229 template ITPP_EXPORT smat transpose(const smat &m);
230 template ITPP_EXPORT imat transpose(const imat &m);
231 template ITPP_EXPORT bmat transpose(const bmat &m);
232 
233 template ITPP_EXPORT void hermitian_transpose(const mat &m, mat &out);
234 template ITPP_EXPORT void hermitian_transpose(const cmat &m, cmat &out);
235 template ITPP_EXPORT void hermitian_transpose(const smat &m, smat &out);
236 template ITPP_EXPORT void hermitian_transpose(const imat &m, imat &out);
237 template ITPP_EXPORT void hermitian_transpose(const bmat &m, bmat &out);
238 
239 template ITPP_EXPORT mat hermitian_transpose(const mat &m);
240 template ITPP_EXPORT cmat hermitian_transpose(const cmat &m);
241 template ITPP_EXPORT smat hermitian_transpose(const smat &m);
242 template ITPP_EXPORT imat hermitian_transpose(const imat &m);
243 template ITPP_EXPORT bmat hermitian_transpose(const bmat &m);
244 
245 template ITPP_EXPORT bool is_hermitian(const mat &X);
246 template ITPP_EXPORT bool is_hermitian(const cmat &X);
247 
248 template ITPP_EXPORT bool is_unitary(const mat &X);
249 template ITPP_EXPORT bool is_unitary(const cmat &X);
250 
251 template ITPP_EXPORT vec rvectorize(const mat &m);
252 template ITPP_EXPORT cvec rvectorize(const cmat &m);
253 template ITPP_EXPORT ivec rvectorize(const imat &m);
254 template ITPP_EXPORT svec rvectorize(const smat &m);
255 template ITPP_EXPORT bvec rvectorize(const bmat &m);
256 
257 template ITPP_EXPORT vec cvectorize(const mat &m);
258 template ITPP_EXPORT cvec cvectorize(const cmat &m);
259 template ITPP_EXPORT ivec cvectorize(const imat &m);
260 template ITPP_EXPORT svec cvectorize(const smat &m);
261 template ITPP_EXPORT bvec cvectorize(const bmat &m);
262 
263 template ITPP_EXPORT mat reshape(const mat &m, int rows, int cols);
264 template ITPP_EXPORT cmat reshape(const cmat &m, int rows, int cols);
265 template ITPP_EXPORT imat reshape(const imat &m, int rows, int cols);
266 template ITPP_EXPORT smat reshape(const smat &m, int rows, int cols);
267 template ITPP_EXPORT bmat reshape(const bmat &m, int rows, int cols);
268 
269 template ITPP_EXPORT mat reshape(const vec &m, int rows, int cols);
270 template ITPP_EXPORT cmat reshape(const cvec &m, int rows, int cols);
271 template ITPP_EXPORT imat reshape(const ivec &m, int rows, int cols);
272 template ITPP_EXPORT smat reshape(const svec &m, int rows, int cols);
273 template ITPP_EXPORT bmat reshape(const bvec &m, int rows, int cols);
274 
275 template ITPP_EXPORT mat kron(const mat &X, const mat &Y);
276 template ITPP_EXPORT cmat kron(const cmat &X, const cmat &Y);
277 template ITPP_EXPORT imat kron(const imat &X, const imat &Y);
278 template ITPP_EXPORT smat kron(const smat &X, const smat &Y);
279 template ITPP_EXPORT bmat kron(const bmat &X, const bmat &Y);
280 
281 template ITPP_EXPORT vec repmat(const vec &v, int n);
282 template ITPP_EXPORT cvec repmat(const cvec &v, int n);
283 template ITPP_EXPORT ivec repmat(const ivec &v, int n);
284 template ITPP_EXPORT svec repmat(const svec &v, int n);
285 template ITPP_EXPORT bvec repmat(const bvec &v, int n);
286 
287 template ITPP_EXPORT mat repmat(const vec &v, int m, int n, bool transpose);
288 template ITPP_EXPORT cmat repmat(const cvec &v, int m, int n, bool transpose);
289 template ITPP_EXPORT imat repmat(const ivec &v, int m, int n, bool transpose);
290 template ITPP_EXPORT smat repmat(const svec &v, int m, int n, bool transpose);
291 template ITPP_EXPORT bmat repmat(const bvec &v, int m, int n, bool transpose);
292 
293 template ITPP_EXPORT mat repmat(const mat &data, int m, int n);
294 template ITPP_EXPORT cmat repmat(const cmat &data, int m, int n);
295 template ITPP_EXPORT imat repmat(const imat &data, int m, int n);
296 template ITPP_EXPORT smat repmat(const smat &data, int m, int n);
297 template ITPP_EXPORT bmat repmat(const bmat &data, int m, int n);
298 
299 } // namespace itpp
300 
SourceForge Logo

Generated on Sat May 25 2013 16:32:19 for IT++ by Doxygen 1.8.2