00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022 template<typename T1>
00023 inline
00024 typename T1::elem_type
00025 accu(const Base<typename T1::elem_type,T1>& X)
00026 {
00027 arma_extra_debug_sigprint();
00028
00029 typedef typename T1::elem_type eT;
00030
00031 const unwrap<T1> tmp(X.get_ref());
00032 const Mat<eT>& A = tmp.M;
00033
00034 const u32 A_n_elem = A.n_elem;
00035 const eT* A_mem = A.mem;
00036
00037 eT val = eT(0);
00038
00039 for(u32 i=0; i<A_n_elem; ++i)
00040 {
00041 val += A_mem[i];
00042 }
00043
00044 return val;
00045 }
00046
00047
00048
00049
00050 template<typename T1>
00051 inline
00052 typename T1::elem_type
00053 accu(const BaseCube<typename T1::elem_type,T1>& X)
00054 {
00055 arma_extra_debug_sigprint();
00056
00057 typedef typename T1::elem_type eT;
00058
00059 const unwrap_cube<T1> tmp(X.get_ref());
00060 const Cube<eT>& A = tmp.M;
00061
00062 const u32 A_n_elem = A.n_elem;
00063 const eT* A_mem = A.mem;
00064
00065 eT val = eT(0);
00066
00067 for(u32 i=0; i<A_n_elem; ++i)
00068 {
00069 val += A_mem[i];
00070 }
00071
00072 return val;
00073 }
00074
00075
00076
00077 #if defined(ARMA_GOOD_COMPILER)
00078
00079
00080
00081
00082 template<typename T1>
00083 inline
00084 typename T1::elem_type
00085 accu(const Op<T1, op_square>& in)
00086 {
00087 arma_extra_debug_sigprint();
00088
00089 typedef typename T1::elem_type eT;
00090
00091 const unwrap<T1> tmp(in.m);
00092 const Mat<eT>& A = tmp.M;
00093
00094 const u32 A_n_elem = A.n_elem;
00095 const eT* A_mem = A.mem;
00096
00097 eT acc = eT(0);
00098
00099 for(u32 i=0; i<A_n_elem; ++i)
00100 {
00101 const eT val = A_mem[i];
00102 acc += val*val;
00103 }
00104
00105 return acc;
00106 }
00107
00108
00109
00110
00111 template<typename T1>
00112 inline
00113 typename T1::elem_type
00114 accu(const OpCube<T1, op_square>& in)
00115 {
00116 arma_extra_debug_sigprint();
00117
00118 typedef typename T1::elem_type eT;
00119
00120 const unwrap_cube<T1> tmp(in.m);
00121 const Cube<eT>& A = tmp.M;
00122
00123 const u32 A_n_elem = A.n_elem;
00124 const eT* A_mem = A.mem;
00125
00126 eT acc = eT(0);
00127
00128 for(u32 i=0; i<A_n_elem; ++i)
00129 {
00130 const eT val = A_mem[i];
00131 acc += val*val;
00132 }
00133
00134 return acc;
00135 }
00136
00137
00138
00139
00140 template<typename T1>
00141 inline
00142 typename T1::elem_type
00143 accu(const Op<T1, op_sqrt>& in)
00144 {
00145 arma_extra_debug_sigprint();
00146
00147 typedef typename T1::elem_type eT;
00148
00149 const unwrap<T1> tmp(in.m);
00150 const Mat<eT>& A = tmp.M;
00151
00152 const u32 A_n_elem = A.n_elem;
00153 const eT* A_mem = A.mem;
00154
00155 eT acc = eT(0);
00156 for(u32 i=0; i<A_n_elem; ++i)
00157 {
00158 acc += std::sqrt(A_mem[i]);
00159 }
00160
00161 return acc;
00162 }
00163
00164
00165
00166
00167 template<typename T1>
00168 inline
00169 typename T1::elem_type
00170 accu(const OpCube<T1, op_sqrt>& in)
00171 {
00172 arma_extra_debug_sigprint();
00173
00174 typedef typename T1::elem_type eT;
00175
00176 const unwrap_cube<T1> tmp(in.m);
00177 const Cube<eT>& A = tmp.M;
00178
00179 const u32 A_n_elem = A.n_elem;
00180 const eT* A_mem = A.mem;
00181
00182 eT acc = eT(0);
00183 for(u32 i=0; i<A_n_elem; ++i)
00184 {
00185 acc += std::sqrt(A_mem[i]);
00186 }
00187
00188 return acc;
00189 }
00190
00191
00192
00193
00194 template<typename T1, typename T2>
00195 inline
00196 typename T1::elem_type
00197 accu(const Op< Glue<T1,T2, glue_minus>, op_square>& in)
00198 {
00199 arma_extra_debug_sigprint();
00200
00201 typedef typename T1::elem_type eT;
00202
00203 const unwrap<T1> tmp1(in.m.A);
00204 const unwrap<T2> tmp2(in.m.B);
00205
00206 const Mat<eT>& A = tmp1.M;
00207 const Mat<eT>& B = tmp2.M;
00208
00209 arma_debug_assert_same_size(A,B, "accu()");
00210
00211 const u32 n_elem = A.n_elem;
00212
00213 eT acc = eT(0);
00214 for(u32 i=0; i<n_elem; ++i)
00215 {
00216 const eT val = A.mem[i] - B.mem[i];
00217 acc += val*val;
00218 }
00219
00220 return acc;
00221 }
00222
00223
00224
00225
00226 template<typename T1, typename T2>
00227 inline
00228 typename T1::elem_type
00229 accu(const OpCube< GlueCube<T1,T2, glue_cube_minus>, op_square>& in)
00230 {
00231 arma_extra_debug_sigprint();
00232
00233 typedef typename T1::elem_type eT;
00234
00235 const unwrap_cube<T1> tmp1(in.m.A);
00236 const unwrap_cube<T2> tmp2(in.m.B);
00237
00238 const Cube<eT>& A = tmp1.M;
00239 const Cube<eT>& B = tmp2.M;
00240
00241 arma_debug_assert_same_size(A,B, "accu()");
00242
00243 const u32 n_elem = A.n_elem;
00244
00245 eT acc = eT(0);
00246 for(u32 i=0; i<n_elem; ++i)
00247 {
00248 const eT val = A.mem[i] - B.mem[i];
00249 acc += val*val;
00250 }
00251
00252 return acc;
00253 }
00254
00255
00256
00257
00258 template<typename eT>
00259 inline
00260 eT
00261 accu_schur(const Mat<eT>& A, const Mat<eT>& B)
00262 {
00263 arma_extra_debug_sigprint();
00264
00265 arma_debug_assert_same_size(A,B, "accu()");
00266
00267 const eT* const A_mem = A.mem;
00268 const eT* const B_mem = B.mem;
00269
00270 const u32 n_elem = A.n_elem;
00271 eT val = eT(0);
00272
00273 for(u32 i=0; i<n_elem; ++i)
00274 {
00275 val += A_mem[i] * B_mem[i];
00276 }
00277
00278 return val;
00279 }
00280
00281
00282
00283
00284 template<typename eT>
00285 inline
00286 eT
00287 accu_schur(const Cube<eT>& A, const Cube<eT>& B)
00288 {
00289 arma_extra_debug_sigprint();
00290
00291 arma_debug_assert_same_size(A,B, "accu()");
00292
00293 const eT* const A_mem = A.mem;
00294 const eT* const B_mem = B.mem;
00295
00296 const u32 n_elem = A.n_elem;
00297 eT val = eT(0);
00298
00299 for(u32 i=0; i<n_elem; ++i)
00300 {
00301 val += A_mem[i] * B_mem[i];
00302 }
00303
00304 return val;
00305 }
00306
00307
00308
00309
00310 template<typename eT>
00311 inline
00312 eT
00313 accu(const Glue<Mat<eT>,Mat<eT>,glue_schur>& X)
00314 {
00315 return accu_schur(X.A, X.B);
00316 }
00317
00318
00319
00320
00321 template<typename eT>
00322 inline
00323 eT
00324 accu(const GlueCube< Cube<eT>, Cube<eT>, glue_cube_schur >& X)
00325 {
00326 return accu_schur(X.A, X.B);
00327 }
00328
00329
00330
00331
00332 template<typename eT>
00333 inline
00334 eT
00335 accu(const Glue<Glue<Mat<eT>,Mat<eT>,glue_schur>,Mat<eT>,glue_schur>& X)
00336 {
00337 arma_extra_debug_sigprint();
00338
00339 const Mat<eT>& A = X.A.A;
00340 const Mat<eT>& B = X.A.B;
00341 const Mat<eT>& C = X.B;
00342
00343 arma_debug_assert_same_size(A,B, "accu()");
00344 arma_debug_assert_same_size(A,C, "accu()");
00345
00346 const eT* const A_mem = A.mem;
00347 const eT* const B_mem = B.mem;
00348 const eT* const C_mem = C.mem;
00349
00350 const u32 n_elem = A.n_elem;
00351 eT val = eT(0);
00352
00353 for(u32 i=0; i<n_elem; ++i)
00354 {
00355 val += A_mem[i] * B_mem[i] * C_mem[i];
00356 }
00357
00358 return val;
00359 }
00360
00361
00362
00363
00364 template<typename eT>
00365 inline
00366 eT
00367 accu(const GlueCube< GlueCube< Cube<eT>, Cube<eT>, glue_cube_schur >, Cube<eT>, glue_cube_schur >& X)
00368 {
00369 arma_extra_debug_sigprint();
00370
00371 const Cube<eT>& A = X.A.A;
00372 const Cube<eT>& B = X.A.B;
00373 const Cube<eT>& C = X.B;
00374
00375 arma_debug_assert_same_size(A,B, "accu()");
00376 arma_debug_assert_same_size(B,C, "accu()");
00377
00378 const eT* const A_mem = A.mem;
00379 const eT* const B_mem = B.mem;
00380 const eT* const C_mem = C.mem;
00381
00382 const u32 n_elem = A.n_elem;
00383 eT val = eT(0);
00384
00385 for(u32 i=0; i<n_elem; ++i)
00386 {
00387 val += A_mem[i] * B_mem[i] * C_mem[i];
00388 }
00389
00390 return val;
00391 }
00392
00393
00394
00395
00396
00397
00398
00399
00400 template<typename T1, typename T2>
00401 inline
00402 typename T1::elem_type
00403 accu(const Glue<T1,T2,glue_schur>& X)
00404 {
00405 arma_extra_debug_sigprint();
00406
00407 isnt_same_type<typename T1::elem_type, typename T2::elem_type>::check();
00408
00409 typedef typename T1::elem_type eT;
00410
00411 const u32 N_mat = 1 + depth_lhs< glue_schur, Glue<T1,T2,glue_schur> >::num;
00412 arma_extra_debug_print( arma_boost::format("N_mat = %d") % N_mat );
00413
00414 if(N_mat == 2)
00415 {
00416 const unwrap<T1> tmp1(X.A);
00417 const unwrap<T2> tmp2(X.B);
00418
00419 return accu_schur(tmp1.M, tmp2.M);
00420 }
00421 else
00422 {
00423 const Mat<eT>* ptrs[N_mat];
00424 bool del[N_mat];
00425
00426 mat_ptrs<glue_schur, Glue<T1,T2,glue_schur> >::get_ptrs(ptrs, del, X);
00427
00428 for(u32 i=0; i<N_mat; ++i) arma_extra_debug_print( arma_boost::format("ptrs[%d] = %x") % i % ptrs[i] );
00429 for(u32 i=0; i<N_mat; ++i) arma_extra_debug_print( arma_boost::format(" del[%d] = %d") % i % del[i] );
00430
00431 const Mat<eT>& tmp_mat = *(ptrs[0]);
00432
00433 for(u32 i=1; i<N_mat; ++i)
00434 {
00435 arma_debug_assert_same_size(tmp_mat, *(ptrs[i]), "accu()");
00436 }
00437
00438
00439
00440
00441 eT val = eT(0);
00442
00443 const u32 n_elem = ptrs[0]->n_elem;
00444
00445 for(u32 j=0; j<n_elem; ++j)
00446 {
00447 eT tmp = ptrs[0]->mem[j];
00448
00449 for(u32 i=1; i<N_mat; ++i)
00450 {
00451 tmp *= ptrs[i]->mem[j];
00452 }
00453
00454 val += tmp;
00455 }
00456
00457
00458 for(u32 i=0; i<N_mat; ++i)
00459 {
00460 if(del[i] == true)
00461 {
00462 arma_extra_debug_print( arma_boost::format("delete mat_ptr[%d]") % i );
00463 delete ptrs[i];
00464 }
00465 }
00466
00467 return val;
00468 }
00469 }
00470
00471
00472
00473
00474 template<typename T1>
00475 inline
00476 typename T1::elem_type
00477 accu(const Op<T1, op_diagmat>& X)
00478 {
00479 arma_extra_debug_sigprint();
00480
00481 typedef typename T1::elem_type eT;
00482
00483 const unwrap<T1> tmp(X.m);
00484 const Mat<eT>& A = tmp.M;
00485
00486 arma_debug_check( !A.is_square(), "accu(): sum of diagonal values of a non-square matrix requested" );
00487
00488 eT acc = eT(0);
00489
00490 for(u32 i=0; i<A.n_rows; ++i)
00491 {
00492 acc += A.at(i,i);
00493 }
00494
00495 return acc;
00496 }
00497
00498
00499
00500 template<typename eT>
00501 inline
00502 eT
00503 accu(const Op<Mat<eT>, op_diagmat_vec>& X)
00504 {
00505 arma_extra_debug_sigprint();
00506
00507 const Mat<eT>& A = X.m;
00508 arma_debug_check( !A.is_vec(), "accu(): internal error: expected a vector" );
00509
00510 return accu(A);
00511 }
00512
00513
00514
00515
00516 template<typename eT>
00517 inline
00518 eT
00519 accu(const diagview<eT>& X)
00520 {
00521 arma_extra_debug_sigprint();
00522
00523 const u32 n_elem = X.n_elem;
00524 eT val = eT(0);
00525
00526 for(u32 i=0; i<n_elem; ++i)
00527 {
00528 val += X[i];
00529 }
00530
00531 return val;
00532 }
00533
00534
00535
00536
00537 template<typename eT>
00538 inline
00539 eT
00540 accu(const subview<eT>& S)
00541 {
00542 arma_extra_debug_sigprint();
00543
00544 eT val = eT(0);
00545
00546 for(u32 col=0; col<S.n_cols; ++col)
00547 {
00548 const eT* coldata = S.colptr(col);
00549
00550 for(u32 row=0; row<S.n_rows; ++row)
00551 {
00552 val += coldata[row];
00553 }
00554
00555 }
00556
00557 return val;
00558 }
00559
00560
00561
00562
00563 template<typename eT>
00564 inline
00565 eT
00566 accu(const subview_row<eT>& S)
00567 {
00568 arma_extra_debug_sigprint();
00569
00570 const Mat<eT>& X = S.m;
00571
00572 const u32 row = S.aux_row1;
00573 const u32 start_col = S.aux_col1;
00574 const u32 end_col = S.aux_col2;
00575
00576 eT val = eT(0);
00577
00578 for(u32 col=start_col; col<=end_col; ++col)
00579 {
00580 val += X.at(row,col);
00581 }
00582
00583 return val;
00584 }
00585
00586
00587
00588
00589 template<typename eT>
00590 inline
00591 eT
00592 accu(const subview_col<eT>& S)
00593 {
00594 arma_extra_debug_sigprint();
00595
00596 const eT* S_colptr = S.colptr(0);
00597 const u32 n_rows = S.n_rows;
00598
00599 eT val = eT(0);
00600
00601 for(u32 row=0; row<n_rows; ++row)
00602 {
00603 val += S_colptr[row];
00604 }
00605
00606 return val;
00607 }
00608
00609
00610
00611
00612 template<typename eT>
00613 inline
00614 eT
00615 accu(const Glue<subview<eT>,Mat<eT>,glue_schur>& X)
00616 {
00617 arma_extra_debug_sigprint();
00618
00619 arma_debug_assert_same_size(X.A, X.B, "accu()");
00620
00621 const Mat<eT>& A = X.A.m;
00622 const Mat<eT>& B = X.B;
00623
00624 const u32 A_sub_n_rows = X.A.n_rows;
00625 const u32 A_sub_n_cols = X.A.n_cols;
00626
00627 const u32 A_aux_row1 = X.A.aux_row1;
00628 const u32 A_aux_col1 = X.A.aux_col1;
00629
00630
00631 eT val = eT(0);
00632
00633 for(u32 col = 0; col<A_sub_n_cols; ++col)
00634 {
00635 const u32 col_mod = A_aux_col1 + col;
00636
00637 for(u32 row = 0; row<A_sub_n_rows; ++row)
00638 {
00639 const u32 row_mod = A_aux_row1 + row;
00640
00641 val += A.at(row_mod, col_mod) * B.at(row,col);
00642 }
00643
00644 }
00645
00646 return val;
00647 }
00648
00649
00650
00651
00652 template<typename eT>
00653 inline
00654 eT
00655 accu(const Glue<Mat<eT>,subview<eT>,glue_schur>& X)
00656 {
00657 arma_extra_debug_sigprint();
00658
00659 arma_debug_assert_same_size(X.A, X.B, "accu()");
00660
00661 const Mat<eT>& A = X.A;
00662 const Mat<eT>& B = X.B.m;
00663
00664
00665
00666
00667 const u32 B_aux_row1 = X.B.aux_row1;
00668 const u32 B_aux_col1 = X.B.aux_col1;
00669
00670
00671 eT val = eT(0);
00672
00673 for(u32 col = 0; col<A.n_cols; ++col)
00674 {
00675 const u32 col_mod = B_aux_col1 + col;
00676
00677 for(u32 row = 0; row<A.n_rows; ++row)
00678 {
00679 const u32 row_mod = B_aux_row1 + row;
00680
00681 val += A.at(row, col) * B.at(row_mod, col_mod);
00682 }
00683
00684 }
00685
00686 return val;
00687 }
00688
00689
00690
00691
00692 template<typename eT>
00693 inline
00694 eT
00695 accu(const Glue<subview<eT>,subview<eT>,glue_schur>& X)
00696 {
00697 arma_extra_debug_sigprint();
00698
00699 arma_debug_assert_same_size(X.A, X.B, "accu()");
00700
00701 const Mat<eT>& A = X.A.m;
00702 const Mat<eT>& B = X.B.m;
00703
00704 const u32 A_sub_n_rows = X.A.n_rows;
00705 const u32 A_sub_n_cols = X.A.n_cols;
00706
00707
00708
00709
00710 const u32 A_aux_row1 = X.A.aux_row1;
00711 const u32 A_aux_col1 = X.A.aux_col1;
00712
00713 const u32 B_aux_row1 = X.B.aux_row1;
00714 const u32 B_aux_col1 = X.B.aux_col1;
00715
00716
00717 eT val = eT(0);
00718
00719 for(u32 col = 0; col<A_sub_n_cols; ++col)
00720 {
00721 const u32 A_col_mod = A_aux_col1 + col;
00722 const u32 B_col_mod = B_aux_col1 + col;
00723
00724 for(u32 row = 0; row<A_sub_n_rows; ++row)
00725 {
00726 const u32 A_row_mod = A_aux_row1 + row;
00727 const u32 B_row_mod = B_aux_row1 + row;
00728
00729 val += A.at(A_row_mod, A_col_mod) * B.at(B_row_mod, B_col_mod);
00730 }
00731
00732 }
00733
00734 return val;
00735 }
00736
00737
00738
00739 #endif
00740
00741
00742
00743