annotate armadillo-3.900.4/include/armadillo_bits/op_htrans_meat.hpp @ 84:55a047986812 tip

Update library URI so as not to be document-local
author Chris Cannam
date Wed, 22 Apr 2020 14:21:57 +0100
parents 1ec0e2823891
children
rev   line source
Chris@49 1 // Copyright (C) 2008-2012 NICTA (www.nicta.com.au)
Chris@49 2 // Copyright (C) 2008-2012 Conrad Sanderson
Chris@49 3 // Copyright (C) 2012 Ryan Curtin
Chris@49 4 //
Chris@49 5 // This Source Code Form is subject to the terms of the Mozilla Public
Chris@49 6 // License, v. 2.0. If a copy of the MPL was not distributed with this
Chris@49 7 // file, You can obtain one at http://mozilla.org/MPL/2.0/.
Chris@49 8
Chris@49 9
Chris@49 10 //! \addtogroup op_htrans
Chris@49 11 //! @{
Chris@49 12
Chris@49 13
Chris@49 14
Chris@49 15 template<typename eT>
Chris@49 16 arma_hot
Chris@49 17 arma_inline
Chris@49 18 void
Chris@49 19 op_htrans::apply_noalias(Mat<eT>& out, const Mat<eT>& A, const typename arma_not_cx<eT>::result* junk)
Chris@49 20 {
Chris@49 21 arma_extra_debug_sigprint();
Chris@49 22 arma_ignore(junk);
Chris@49 23
Chris@49 24 op_strans::apply_noalias(out, A);
Chris@49 25 }
Chris@49 26
Chris@49 27
Chris@49 28
Chris@49 29 template<typename eT>
Chris@49 30 arma_hot
Chris@49 31 inline
Chris@49 32 void
Chris@49 33 op_htrans::apply_noalias(Mat<eT>& out, const Mat<eT>& A, const typename arma_cx_only<eT>::result* junk)
Chris@49 34 {
Chris@49 35 arma_extra_debug_sigprint();
Chris@49 36 arma_ignore(junk);
Chris@49 37
Chris@49 38 const uword A_n_rows = A.n_rows;
Chris@49 39 const uword A_n_cols = A.n_cols;
Chris@49 40
Chris@49 41 out.set_size(A_n_cols, A_n_rows);
Chris@49 42
Chris@49 43 if( (A_n_cols == 1) || (A_n_rows == 1) )
Chris@49 44 {
Chris@49 45 const uword n_elem = A.n_elem;
Chris@49 46
Chris@49 47 const eT* A_mem = A.memptr();
Chris@49 48 eT* out_mem = out.memptr();
Chris@49 49
Chris@49 50 for(uword i=0; i < n_elem; ++i)
Chris@49 51 {
Chris@49 52 out_mem[i] = std::conj(A_mem[i]);
Chris@49 53 }
Chris@49 54 }
Chris@49 55 else
Chris@49 56 {
Chris@49 57 for(uword in_row = 0; in_row < A_n_rows; ++in_row)
Chris@49 58 {
Chris@49 59 const uword out_col = in_row;
Chris@49 60
Chris@49 61 for(uword in_col = 0; in_col < A_n_cols; ++in_col)
Chris@49 62 {
Chris@49 63 const uword out_row = in_col;
Chris@49 64 out.at(out_row, out_col) = std::conj( A.at(in_row, in_col) );
Chris@49 65 }
Chris@49 66 }
Chris@49 67 }
Chris@49 68 }
Chris@49 69
Chris@49 70
Chris@49 71
Chris@49 72 template<typename eT>
Chris@49 73 arma_hot
Chris@49 74 arma_inline
Chris@49 75 void
Chris@49 76 op_htrans::apply(Mat<eT>& out, const Mat<eT>& A, const typename arma_not_cx<eT>::result* junk)
Chris@49 77 {
Chris@49 78 arma_extra_debug_sigprint();
Chris@49 79 arma_ignore(junk);
Chris@49 80
Chris@49 81 op_strans::apply(out, A);
Chris@49 82 }
Chris@49 83
Chris@49 84
Chris@49 85
Chris@49 86 template<typename eT>
Chris@49 87 arma_hot
Chris@49 88 inline
Chris@49 89 void
Chris@49 90 op_htrans::apply(Mat<eT>& out, const Mat<eT>& A, const typename arma_cx_only<eT>::result* junk)
Chris@49 91 {
Chris@49 92 arma_extra_debug_sigprint();
Chris@49 93 arma_ignore(junk);
Chris@49 94
Chris@49 95 if(&out != &A)
Chris@49 96 {
Chris@49 97 op_htrans::apply_noalias(out, A);
Chris@49 98 }
Chris@49 99 else
Chris@49 100 {
Chris@49 101 const uword n_rows = out.n_rows;
Chris@49 102 const uword n_cols = out.n_cols;
Chris@49 103
Chris@49 104 if(n_rows == n_cols)
Chris@49 105 {
Chris@49 106 arma_extra_debug_print("doing in-place hermitian transpose of a square matrix");
Chris@49 107
Chris@49 108 for(uword col=0; col < n_cols; ++col)
Chris@49 109 {
Chris@49 110 eT* coldata = out.colptr(col);
Chris@49 111
Chris@49 112 out.at(col,col) = std::conj( out.at(col,col) );
Chris@49 113
Chris@49 114 for(uword row=(col+1); row < n_rows; ++row)
Chris@49 115 {
Chris@49 116 const eT val1 = std::conj(coldata[row]);
Chris@49 117 const eT val2 = std::conj(out.at(col,row));
Chris@49 118
Chris@49 119 out.at(col,row) = val1;
Chris@49 120 coldata[row] = val2;
Chris@49 121 }
Chris@49 122 }
Chris@49 123 }
Chris@49 124 else
Chris@49 125 {
Chris@49 126 Mat<eT> tmp;
Chris@49 127 op_htrans::apply_noalias(tmp, A);
Chris@49 128
Chris@49 129 out.steal_mem(tmp);
Chris@49 130 }
Chris@49 131 }
Chris@49 132
Chris@49 133 }
Chris@49 134
Chris@49 135
Chris@49 136
Chris@49 137 template<typename T1>
Chris@49 138 arma_hot
Chris@49 139 inline
Chris@49 140 void
Chris@49 141 op_htrans::apply_proxy(Mat<typename T1::elem_type>& out, const T1& X)
Chris@49 142 {
Chris@49 143 arma_extra_debug_sigprint();
Chris@49 144
Chris@49 145 typedef typename T1::elem_type eT;
Chris@49 146
Chris@49 147 const Proxy<T1> P(X);
Chris@49 148
Chris@49 149 // allow detection of in-place transpose
Chris@49 150 if( (is_Mat<typename Proxy<T1>::stored_type>::value == true) && (Proxy<T1>::fake_mat == false) )
Chris@49 151 {
Chris@49 152 const unwrap<typename Proxy<T1>::stored_type> tmp(P.Q);
Chris@49 153
Chris@49 154 op_htrans::apply(out, tmp.M);
Chris@49 155 }
Chris@49 156 else
Chris@49 157 {
Chris@49 158 const uword n_rows = P.get_n_rows();
Chris@49 159 const uword n_cols = P.get_n_cols();
Chris@49 160
Chris@49 161 const bool is_alias = P.is_alias(out);
Chris@49 162
Chris@49 163 if( (resolves_to_vector<T1>::value == true) && (Proxy<T1>::prefer_at_accessor == false) )
Chris@49 164 {
Chris@49 165 if(is_alias == false)
Chris@49 166 {
Chris@49 167 out.set_size(n_cols, n_rows);
Chris@49 168
Chris@49 169 eT* out_mem = out.memptr();
Chris@49 170
Chris@49 171 const uword n_elem = P.get_n_elem();
Chris@49 172
Chris@49 173 typename Proxy<T1>::ea_type Pea = P.get_ea();
Chris@49 174
Chris@49 175 for(uword i=0; i < n_elem; ++i)
Chris@49 176 {
Chris@49 177 out_mem[i] = std::conj(Pea[i]);
Chris@49 178 }
Chris@49 179 }
Chris@49 180 else // aliasing
Chris@49 181 {
Chris@49 182 Mat<eT> out2(n_cols, n_rows);
Chris@49 183
Chris@49 184 eT* out_mem = out2.memptr();
Chris@49 185
Chris@49 186 const uword n_elem = P.get_n_elem();
Chris@49 187
Chris@49 188 typename Proxy<T1>::ea_type Pea = P.get_ea();
Chris@49 189
Chris@49 190 for(uword i=0; i < n_elem; ++i)
Chris@49 191 {
Chris@49 192 out_mem[i] = std::conj(Pea[i]);
Chris@49 193 }
Chris@49 194
Chris@49 195 out.steal_mem(out2);
Chris@49 196 }
Chris@49 197 }
Chris@49 198 else
Chris@49 199 {
Chris@49 200 if(is_alias == false)
Chris@49 201 {
Chris@49 202 out.set_size(n_cols, n_rows);
Chris@49 203
Chris@49 204 for(uword k=0; k < n_cols; ++k)
Chris@49 205 for(uword i=0; i < n_rows; ++i)
Chris@49 206 {
Chris@49 207 out.at(k,i) = std::conj(P.at(i,k));
Chris@49 208 }
Chris@49 209 }
Chris@49 210 else // aliasing
Chris@49 211 {
Chris@49 212 Mat<eT> out2(n_cols, n_rows);
Chris@49 213
Chris@49 214 for(uword k=0; k < n_cols; ++k)
Chris@49 215 for(uword i=0; i < n_rows; ++i)
Chris@49 216 {
Chris@49 217 out2.at(k,i) = std::conj(P.at(i,k));
Chris@49 218 }
Chris@49 219
Chris@49 220 out.steal_mem(out2);
Chris@49 221 }
Chris@49 222 }
Chris@49 223 }
Chris@49 224 }
Chris@49 225
Chris@49 226
Chris@49 227
Chris@49 228 template<typename T1>
Chris@49 229 arma_hot
Chris@49 230 inline
Chris@49 231 void
Chris@49 232 op_htrans::apply(Mat<typename T1::elem_type>& out, const Op<T1,op_htrans>& in, const typename arma_not_cx<typename T1::elem_type>::result* junk)
Chris@49 233 {
Chris@49 234 arma_extra_debug_sigprint();
Chris@49 235 arma_ignore(junk);
Chris@49 236
Chris@49 237 op_strans::apply_proxy(out, in.m);
Chris@49 238 }
Chris@49 239
Chris@49 240
Chris@49 241
Chris@49 242 template<typename T1>
Chris@49 243 arma_hot
Chris@49 244 inline
Chris@49 245 void
Chris@49 246 op_htrans::apply(Mat<typename T1::elem_type>& out, const Op<T1,op_htrans>& in, const typename arma_cx_only<typename T1::elem_type>::result* junk)
Chris@49 247 {
Chris@49 248 arma_extra_debug_sigprint();
Chris@49 249 arma_ignore(junk);
Chris@49 250
Chris@49 251 op_htrans::apply_proxy(out, in.m);
Chris@49 252 }
Chris@49 253
Chris@49 254
Chris@49 255
Chris@49 256 template<typename T1>
Chris@49 257 arma_hot
Chris@49 258 inline
Chris@49 259 void
Chris@49 260 op_htrans::apply(Mat<typename T1::elem_type>& out, const Op< Op<T1, op_trimat>, op_htrans>& in)
Chris@49 261 {
Chris@49 262 arma_extra_debug_sigprint();
Chris@49 263
Chris@49 264 typedef typename T1::elem_type eT;
Chris@49 265
Chris@49 266 const unwrap<T1> tmp(in.m.m);
Chris@49 267 const Mat<eT>& A = tmp.M;
Chris@49 268
Chris@49 269 const bool upper = in.m.aux_uword_a;
Chris@49 270
Chris@49 271 op_trimat::apply_htrans(out, A, upper);
Chris@49 272 }
Chris@49 273
Chris@49 274
Chris@49 275
Chris@49 276 //
Chris@49 277 // op_htrans2
Chris@49 278
Chris@49 279
Chris@49 280
Chris@49 281 template<typename eT>
Chris@49 282 arma_hot
Chris@49 283 arma_inline
Chris@49 284 void
Chris@49 285 op_htrans2::apply_noalias(Mat<eT>& out, const Mat<eT>& A, const eT val)
Chris@49 286 {
Chris@49 287 arma_extra_debug_sigprint();
Chris@49 288
Chris@49 289 const uword A_n_rows = A.n_rows;
Chris@49 290 const uword A_n_cols = A.n_cols;
Chris@49 291
Chris@49 292 out.set_size(A_n_cols, A_n_rows);
Chris@49 293
Chris@49 294 if( (A_n_cols == 1) || (A_n_rows == 1) )
Chris@49 295 {
Chris@49 296 const uword n_elem = A.n_elem;
Chris@49 297
Chris@49 298 const eT* A_mem = A.memptr();
Chris@49 299 eT* out_mem = out.memptr();
Chris@49 300
Chris@49 301 for(uword i=0; i < n_elem; ++i)
Chris@49 302 {
Chris@49 303 out_mem[i] = val * std::conj(A_mem[i]);
Chris@49 304 }
Chris@49 305 }
Chris@49 306 else
Chris@49 307 {
Chris@49 308 for(uword in_row = 0; in_row < A_n_rows; ++in_row)
Chris@49 309 {
Chris@49 310 const uword out_col = in_row;
Chris@49 311
Chris@49 312 for(uword in_col = 0; in_col < A_n_cols; ++in_col)
Chris@49 313 {
Chris@49 314 const uword out_row = in_col;
Chris@49 315 out.at(out_row, out_col) = val * std::conj( A.at(in_row, in_col) );
Chris@49 316 }
Chris@49 317 }
Chris@49 318 }
Chris@49 319 }
Chris@49 320
Chris@49 321
Chris@49 322
Chris@49 323 template<typename eT>
Chris@49 324 arma_hot
Chris@49 325 inline
Chris@49 326 void
Chris@49 327 op_htrans2::apply(Mat<eT>& out, const Mat<eT>& A, const eT val)
Chris@49 328 {
Chris@49 329 arma_extra_debug_sigprint();
Chris@49 330
Chris@49 331 if(&out != &A)
Chris@49 332 {
Chris@49 333 op_htrans2::apply_noalias(out, A, val);
Chris@49 334 }
Chris@49 335 else
Chris@49 336 {
Chris@49 337 const uword n_rows = out.n_rows;
Chris@49 338 const uword n_cols = out.n_cols;
Chris@49 339
Chris@49 340 if(n_rows == n_cols)
Chris@49 341 {
Chris@49 342 arma_extra_debug_print("doing in-place hermitian transpose of a square matrix");
Chris@49 343
Chris@49 344 // TODO: do multiplication while swapping
Chris@49 345
Chris@49 346 for(uword col=0; col < n_cols; ++col)
Chris@49 347 {
Chris@49 348 eT* coldata = out.colptr(col);
Chris@49 349
Chris@49 350 out.at(col,col) = std::conj( out.at(col,col) );
Chris@49 351
Chris@49 352 for(uword row=(col+1); row < n_rows; ++row)
Chris@49 353 {
Chris@49 354 const eT val1 = std::conj(coldata[row]);
Chris@49 355 const eT val2 = std::conj(out.at(col,row));
Chris@49 356
Chris@49 357 out.at(col,row) = val1;
Chris@49 358 coldata[row] = val2;
Chris@49 359 }
Chris@49 360 }
Chris@49 361
Chris@49 362 arrayops::inplace_mul( out.memptr(), val, out.n_elem );
Chris@49 363 }
Chris@49 364 else
Chris@49 365 {
Chris@49 366 Mat<eT> tmp;
Chris@49 367 op_htrans2::apply_noalias(tmp, A, val);
Chris@49 368
Chris@49 369 out.steal_mem(tmp);
Chris@49 370 }
Chris@49 371 }
Chris@49 372 }
Chris@49 373
Chris@49 374
Chris@49 375
Chris@49 376 template<typename T1>
Chris@49 377 arma_hot
Chris@49 378 inline
Chris@49 379 void
Chris@49 380 op_htrans2::apply_proxy(Mat<typename T1::elem_type>& out, const T1& X, const typename T1::elem_type val)
Chris@49 381 {
Chris@49 382 arma_extra_debug_sigprint();
Chris@49 383
Chris@49 384 typedef typename T1::elem_type eT;
Chris@49 385
Chris@49 386 const Proxy<T1> P(X);
Chris@49 387
Chris@49 388 // allow detection of in-place transpose
Chris@49 389 if( (is_Mat<typename Proxy<T1>::stored_type>::value == true) && (Proxy<T1>::fake_mat == false) )
Chris@49 390 {
Chris@49 391 const unwrap<typename Proxy<T1>::stored_type> tmp(P.Q);
Chris@49 392
Chris@49 393 op_htrans2::apply(out, tmp.M, val);
Chris@49 394 }
Chris@49 395 else
Chris@49 396 {
Chris@49 397 const uword n_rows = P.get_n_rows();
Chris@49 398 const uword n_cols = P.get_n_cols();
Chris@49 399
Chris@49 400 const bool is_alias = P.is_alias(out);
Chris@49 401
Chris@49 402 if( (resolves_to_vector<T1>::value == true) && (Proxy<T1>::prefer_at_accessor == false) )
Chris@49 403 {
Chris@49 404 if(is_alias == false)
Chris@49 405 {
Chris@49 406 out.set_size(n_cols, n_rows);
Chris@49 407
Chris@49 408 eT* out_mem = out.memptr();
Chris@49 409
Chris@49 410 const uword n_elem = P.get_n_elem();
Chris@49 411
Chris@49 412 typename Proxy<T1>::ea_type Pea = P.get_ea();
Chris@49 413
Chris@49 414 for(uword i=0; i < n_elem; ++i)
Chris@49 415 {
Chris@49 416 out_mem[i] = val * std::conj(Pea[i]);
Chris@49 417 }
Chris@49 418 }
Chris@49 419 else // aliasing
Chris@49 420 {
Chris@49 421 Mat<eT> out2(n_cols, n_rows);
Chris@49 422
Chris@49 423 eT* out_mem = out2.memptr();
Chris@49 424
Chris@49 425 const uword n_elem = P.get_n_elem();
Chris@49 426
Chris@49 427 typename Proxy<T1>::ea_type Pea = P.get_ea();
Chris@49 428
Chris@49 429 for(uword i=0; i < n_elem; ++i)
Chris@49 430 {
Chris@49 431 out_mem[i] = val * std::conj(Pea[i]);
Chris@49 432 }
Chris@49 433
Chris@49 434 out.steal_mem(out2);
Chris@49 435 }
Chris@49 436 }
Chris@49 437 else
Chris@49 438 {
Chris@49 439 if(is_alias == false)
Chris@49 440 {
Chris@49 441 out.set_size(n_cols, n_rows);
Chris@49 442
Chris@49 443 for(uword k=0; k < n_cols; ++k)
Chris@49 444 for(uword i=0; i < n_rows; ++i)
Chris@49 445 {
Chris@49 446 out.at(k,i) = val * std::conj(P.at(i,k));
Chris@49 447 }
Chris@49 448 }
Chris@49 449 else // aliasing
Chris@49 450 {
Chris@49 451 Mat<eT> out2(n_cols, n_rows);
Chris@49 452
Chris@49 453 for(uword k=0; k < n_cols; ++k)
Chris@49 454 for(uword i=0; i < n_rows; ++i)
Chris@49 455 {
Chris@49 456 out2.at(k,i) = val * std::conj(P.at(i,k));
Chris@49 457 }
Chris@49 458
Chris@49 459 out.steal_mem(out2);
Chris@49 460 }
Chris@49 461 }
Chris@49 462 }
Chris@49 463 }
Chris@49 464
Chris@49 465
Chris@49 466
Chris@49 467 template<typename T1>
Chris@49 468 arma_hot
Chris@49 469 inline
Chris@49 470 void
Chris@49 471 op_htrans2::apply(Mat<typename T1::elem_type>& out, const Op<T1,op_htrans2>& in, const typename arma_not_cx<typename T1::elem_type>::result* junk)
Chris@49 472 {
Chris@49 473 arma_extra_debug_sigprint();
Chris@49 474 arma_ignore(junk);
Chris@49 475
Chris@49 476 op_strans2::apply_proxy(out, in.m, in.aux);
Chris@49 477 }
Chris@49 478
Chris@49 479
Chris@49 480
Chris@49 481 template<typename T1>
Chris@49 482 arma_hot
Chris@49 483 inline
Chris@49 484 void
Chris@49 485 op_htrans2::apply(Mat<typename T1::elem_type>& out, const Op<T1,op_htrans2>& in, const typename arma_cx_only<typename T1::elem_type>::result* junk)
Chris@49 486 {
Chris@49 487 arma_extra_debug_sigprint();
Chris@49 488 arma_ignore(junk);
Chris@49 489
Chris@49 490 op_htrans2::apply_proxy(out, in.m, in.aux);
Chris@49 491 }
Chris@49 492
Chris@49 493
Chris@49 494
Chris@49 495 //! @}