diagview_meat.hpp 21 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840841842843844845846847848849850851852853854855856857858859860861862863864865866867868869870871872873874875876877878879880881882883884885886887888889890891892893894895896897898899900901902903904905906907908909910911912913914915916917918919920921922923924925926927928929930931932933934935936937938939940941942943944945946947948949950951952953954955956957958959960961962963964965966967968969970971972973974975976977978979980981982983984985
  1. // Copyright 2008-2016 Conrad Sanderson (http://conradsanderson.id.au)
  2. // Copyright 2008-2016 National ICT Australia (NICTA)
  3. //
  4. // Licensed under the Apache License, Version 2.0 (the "License");
  5. // you may not use this file except in compliance with the License.
  6. // You may obtain a copy of the License at
  7. // http://www.apache.org/licenses/LICENSE-2.0
  8. //
  9. // Unless required by applicable law or agreed to in writing, software
  10. // distributed under the License is distributed on an "AS IS" BASIS,
  11. // WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
  12. // See the License for the specific language governing permissions and
  13. // limitations under the License.
  14. // ------------------------------------------------------------------------
  15. //! \addtogroup diagview
  16. //! @{
  17. template<typename eT>
  18. inline
  19. diagview<eT>::~diagview()
  20. {
  21. arma_extra_debug_sigprint();
  22. }
  23. template<typename eT>
  24. arma_inline
  25. diagview<eT>::diagview(const Mat<eT>& in_m, const uword in_row_offset, const uword in_col_offset, const uword in_len)
  26. : m(in_m)
  27. , row_offset(in_row_offset)
  28. , col_offset(in_col_offset)
  29. , n_rows(in_len)
  30. , n_elem(in_len)
  31. {
  32. arma_extra_debug_sigprint();
  33. }
  34. //! set a diagonal of our matrix using a diagonal from a foreign matrix
  35. template<typename eT>
  36. inline
  37. void
  38. diagview<eT>::operator= (const diagview<eT>& x)
  39. {
  40. arma_extra_debug_sigprint();
  41. diagview<eT>& d = *this;
  42. arma_debug_check( (d.n_elem != x.n_elem), "diagview: diagonals have incompatible lengths");
  43. Mat<eT>& d_m = const_cast< Mat<eT>& >(d.m);
  44. const Mat<eT>& x_m = x.m;
  45. if(&d_m != &x_m)
  46. {
  47. const uword d_n_elem = d.n_elem;
  48. const uword d_row_offset = d.row_offset;
  49. const uword d_col_offset = d.col_offset;
  50. const uword x_row_offset = x.row_offset;
  51. const uword x_col_offset = x.col_offset;
  52. uword ii,jj;
  53. for(ii=0, jj=1; jj < d_n_elem; ii+=2, jj+=2)
  54. {
  55. const eT tmp_i = x_m.at(ii + x_row_offset, ii + x_col_offset);
  56. const eT tmp_j = x_m.at(jj + x_row_offset, jj + x_col_offset);
  57. d_m.at(ii + d_row_offset, ii + d_col_offset) = tmp_i;
  58. d_m.at(jj + d_row_offset, jj + d_col_offset) = tmp_j;
  59. }
  60. if(ii < d_n_elem)
  61. {
  62. d_m.at(ii + d_row_offset, ii + d_col_offset) = x_m.at(ii + x_row_offset, ii + x_col_offset);
  63. }
  64. }
  65. else
  66. {
  67. const Mat<eT> tmp = x;
  68. (*this).operator=(tmp);
  69. }
  70. }
  71. template<typename eT>
  72. inline
  73. void
  74. diagview<eT>::operator+=(const eT val)
  75. {
  76. arma_extra_debug_sigprint();
  77. Mat<eT>& t_m = const_cast< Mat<eT>& >(m);
  78. const uword t_n_elem = n_elem;
  79. const uword t_row_offset = row_offset;
  80. const uword t_col_offset = col_offset;
  81. for(uword ii=0; ii < t_n_elem; ++ii)
  82. {
  83. t_m.at( ii + t_row_offset, ii + t_col_offset) += val;
  84. }
  85. }
  86. template<typename eT>
  87. inline
  88. void
  89. diagview<eT>::operator-=(const eT val)
  90. {
  91. arma_extra_debug_sigprint();
  92. Mat<eT>& t_m = const_cast< Mat<eT>& >(m);
  93. const uword t_n_elem = n_elem;
  94. const uword t_row_offset = row_offset;
  95. const uword t_col_offset = col_offset;
  96. for(uword ii=0; ii < t_n_elem; ++ii)
  97. {
  98. t_m.at( ii + t_row_offset, ii + t_col_offset) -= val;
  99. }
  100. }
  101. template<typename eT>
  102. inline
  103. void
  104. diagview<eT>::operator*=(const eT val)
  105. {
  106. arma_extra_debug_sigprint();
  107. Mat<eT>& t_m = const_cast< Mat<eT>& >(m);
  108. const uword t_n_elem = n_elem;
  109. const uword t_row_offset = row_offset;
  110. const uword t_col_offset = col_offset;
  111. for(uword ii=0; ii < t_n_elem; ++ii)
  112. {
  113. t_m.at( ii + t_row_offset, ii + t_col_offset) *= val;
  114. }
  115. }
  116. template<typename eT>
  117. inline
  118. void
  119. diagview<eT>::operator/=(const eT val)
  120. {
  121. arma_extra_debug_sigprint();
  122. Mat<eT>& t_m = const_cast< Mat<eT>& >(m);
  123. const uword t_n_elem = n_elem;
  124. const uword t_row_offset = row_offset;
  125. const uword t_col_offset = col_offset;
  126. for(uword ii=0; ii < t_n_elem; ++ii)
  127. {
  128. t_m.at( ii + t_row_offset, ii + t_col_offset) /= val;
  129. }
  130. }
  131. //! set a diagonal of our matrix using data from a foreign object
  132. template<typename eT>
  133. template<typename T1>
  134. inline
  135. void
  136. diagview<eT>::operator= (const Base<eT,T1>& o)
  137. {
  138. arma_extra_debug_sigprint();
  139. diagview<eT>& d = *this;
  140. Mat<eT>& d_m = const_cast< Mat<eT>& >(d.m);
  141. const uword d_n_elem = d.n_elem;
  142. const uword d_row_offset = d.row_offset;
  143. const uword d_col_offset = d.col_offset;
  144. const Proxy<T1> P( o.get_ref() );
  145. arma_debug_check
  146. (
  147. ( (d_n_elem != P.get_n_elem()) || ((P.get_n_rows() != 1) && (P.get_n_cols() != 1)) ),
  148. "diagview: given object has incompatible size"
  149. );
  150. const bool is_alias = P.is_alias(d_m);
  151. if(is_alias) { arma_extra_debug_print("aliasing detected"); }
  152. if( (is_Mat<typename Proxy<T1>::stored_type>::value) || (Proxy<T1>::use_at) || (is_alias) )
  153. {
  154. const unwrap_check<typename Proxy<T1>::stored_type> tmp(P.Q, is_alias);
  155. const Mat<eT>& x = tmp.M;
  156. const eT* x_mem = x.memptr();
  157. uword ii,jj;
  158. for(ii=0, jj=1; jj < d_n_elem; ii+=2, jj+=2)
  159. {
  160. const eT tmp_i = x_mem[ii];
  161. const eT tmp_j = x_mem[jj];
  162. d_m.at( ii + d_row_offset, ii + d_col_offset) = tmp_i;
  163. d_m.at( jj + d_row_offset, jj + d_col_offset) = tmp_j;
  164. }
  165. if(ii < d_n_elem)
  166. {
  167. d_m.at( ii + d_row_offset, ii + d_col_offset) = x_mem[ii];
  168. }
  169. }
  170. else
  171. {
  172. typename Proxy<T1>::ea_type Pea = P.get_ea();
  173. uword ii,jj;
  174. for(ii=0, jj=1; jj < d_n_elem; ii+=2, jj+=2)
  175. {
  176. const eT tmp_i = Pea[ii];
  177. const eT tmp_j = Pea[jj];
  178. d_m.at( ii + d_row_offset, ii + d_col_offset) = tmp_i;
  179. d_m.at( jj + d_row_offset, jj + d_col_offset) = tmp_j;
  180. }
  181. if(ii < d_n_elem)
  182. {
  183. d_m.at( ii + d_row_offset, ii + d_col_offset) = Pea[ii];
  184. }
  185. }
  186. }
  187. template<typename eT>
  188. template<typename T1>
  189. inline
  190. void
  191. diagview<eT>::operator+=(const Base<eT,T1>& o)
  192. {
  193. arma_extra_debug_sigprint();
  194. diagview<eT>& d = *this;
  195. Mat<eT>& d_m = const_cast< Mat<eT>& >(d.m);
  196. const uword d_n_elem = d.n_elem;
  197. const uword d_row_offset = d.row_offset;
  198. const uword d_col_offset = d.col_offset;
  199. const Proxy<T1> P( o.get_ref() );
  200. arma_debug_check
  201. (
  202. ( (d_n_elem != P.get_n_elem()) || ((P.get_n_rows() != 1) && (P.get_n_cols() != 1)) ),
  203. "diagview: given object has incompatible size"
  204. );
  205. const bool is_alias = P.is_alias(d_m);
  206. if(is_alias) { arma_extra_debug_print("aliasing detected"); }
  207. if( (is_Mat<typename Proxy<T1>::stored_type>::value) || (Proxy<T1>::use_at) || (is_alias) )
  208. {
  209. const unwrap_check<typename Proxy<T1>::stored_type> tmp(P.Q, is_alias);
  210. const Mat<eT>& x = tmp.M;
  211. const eT* x_mem = x.memptr();
  212. uword ii,jj;
  213. for(ii=0, jj=1; jj < d_n_elem; ii+=2, jj+=2)
  214. {
  215. const eT tmp_i = x_mem[ii];
  216. const eT tmp_j = x_mem[jj];
  217. d_m.at( ii + d_row_offset, ii + d_col_offset) += tmp_i;
  218. d_m.at( jj + d_row_offset, jj + d_col_offset) += tmp_j;
  219. }
  220. if(ii < d_n_elem)
  221. {
  222. d_m.at( ii + d_row_offset, ii + d_col_offset) += x_mem[ii];
  223. }
  224. }
  225. else
  226. {
  227. typename Proxy<T1>::ea_type Pea = P.get_ea();
  228. uword ii,jj;
  229. for(ii=0, jj=1; jj < d_n_elem; ii+=2, jj+=2)
  230. {
  231. const eT tmp_i = Pea[ii];
  232. const eT tmp_j = Pea[jj];
  233. d_m.at( ii + d_row_offset, ii + d_col_offset) += tmp_i;
  234. d_m.at( jj + d_row_offset, jj + d_col_offset) += tmp_j;
  235. }
  236. if(ii < d_n_elem)
  237. {
  238. d_m.at( ii + d_row_offset, ii + d_col_offset) += Pea[ii];
  239. }
  240. }
  241. }
  242. template<typename eT>
  243. template<typename T1>
  244. inline
  245. void
  246. diagview<eT>::operator-=(const Base<eT,T1>& o)
  247. {
  248. arma_extra_debug_sigprint();
  249. diagview<eT>& d = *this;
  250. Mat<eT>& d_m = const_cast< Mat<eT>& >(d.m);
  251. const uword d_n_elem = d.n_elem;
  252. const uword d_row_offset = d.row_offset;
  253. const uword d_col_offset = d.col_offset;
  254. const Proxy<T1> P( o.get_ref() );
  255. arma_debug_check
  256. (
  257. ( (d_n_elem != P.get_n_elem()) || ((P.get_n_rows() != 1) && (P.get_n_cols() != 1)) ),
  258. "diagview: given object has incompatible size"
  259. );
  260. const bool is_alias = P.is_alias(d_m);
  261. if(is_alias) { arma_extra_debug_print("aliasing detected"); }
  262. if( (is_Mat<typename Proxy<T1>::stored_type>::value) || (Proxy<T1>::use_at) || (is_alias) )
  263. {
  264. const unwrap_check<typename Proxy<T1>::stored_type> tmp(P.Q, is_alias);
  265. const Mat<eT>& x = tmp.M;
  266. const eT* x_mem = x.memptr();
  267. uword ii,jj;
  268. for(ii=0, jj=1; jj < d_n_elem; ii+=2, jj+=2)
  269. {
  270. const eT tmp_i = x_mem[ii];
  271. const eT tmp_j = x_mem[jj];
  272. d_m.at( ii + d_row_offset, ii + d_col_offset) -= tmp_i;
  273. d_m.at( jj + d_row_offset, jj + d_col_offset) -= tmp_j;
  274. }
  275. if(ii < d_n_elem)
  276. {
  277. d_m.at( ii + d_row_offset, ii + d_col_offset) -= x_mem[ii];
  278. }
  279. }
  280. else
  281. {
  282. typename Proxy<T1>::ea_type Pea = P.get_ea();
  283. uword ii,jj;
  284. for(ii=0, jj=1; jj < d_n_elem; ii+=2, jj+=2)
  285. {
  286. const eT tmp_i = Pea[ii];
  287. const eT tmp_j = Pea[jj];
  288. d_m.at( ii + d_row_offset, ii + d_col_offset) -= tmp_i;
  289. d_m.at( jj + d_row_offset, jj + d_col_offset) -= tmp_j;
  290. }
  291. if(ii < d_n_elem)
  292. {
  293. d_m.at( ii + d_row_offset, ii + d_col_offset) -= Pea[ii];
  294. }
  295. }
  296. }
  297. template<typename eT>
  298. template<typename T1>
  299. inline
  300. void
  301. diagview<eT>::operator%=(const Base<eT,T1>& o)
  302. {
  303. arma_extra_debug_sigprint();
  304. diagview<eT>& d = *this;
  305. Mat<eT>& d_m = const_cast< Mat<eT>& >(d.m);
  306. const uword d_n_elem = d.n_elem;
  307. const uword d_row_offset = d.row_offset;
  308. const uword d_col_offset = d.col_offset;
  309. const Proxy<T1> P( o.get_ref() );
  310. arma_debug_check
  311. (
  312. ( (d_n_elem != P.get_n_elem()) || ((P.get_n_rows() != 1) && (P.get_n_cols() != 1)) ),
  313. "diagview: given object has incompatible size"
  314. );
  315. const bool is_alias = P.is_alias(d_m);
  316. if(is_alias) { arma_extra_debug_print("aliasing detected"); }
  317. if( (is_Mat<typename Proxy<T1>::stored_type>::value) || (Proxy<T1>::use_at) || (is_alias) )
  318. {
  319. const unwrap_check<typename Proxy<T1>::stored_type> tmp(P.Q, is_alias);
  320. const Mat<eT>& x = tmp.M;
  321. const eT* x_mem = x.memptr();
  322. uword ii,jj;
  323. for(ii=0, jj=1; jj < d_n_elem; ii+=2, jj+=2)
  324. {
  325. const eT tmp_i = x_mem[ii];
  326. const eT tmp_j = x_mem[jj];
  327. d_m.at( ii + d_row_offset, ii + d_col_offset) *= tmp_i;
  328. d_m.at( jj + d_row_offset, jj + d_col_offset) *= tmp_j;
  329. }
  330. if(ii < d_n_elem)
  331. {
  332. d_m.at( ii + d_row_offset, ii + d_col_offset) *= x_mem[ii];
  333. }
  334. }
  335. else
  336. {
  337. typename Proxy<T1>::ea_type Pea = P.get_ea();
  338. uword ii,jj;
  339. for(ii=0, jj=1; jj < d_n_elem; ii+=2, jj+=2)
  340. {
  341. const eT tmp_i = Pea[ii];
  342. const eT tmp_j = Pea[jj];
  343. d_m.at( ii + d_row_offset, ii + d_col_offset) *= tmp_i;
  344. d_m.at( jj + d_row_offset, jj + d_col_offset) *= tmp_j;
  345. }
  346. if(ii < d_n_elem)
  347. {
  348. d_m.at( ii + d_row_offset, ii + d_col_offset) *= Pea[ii];
  349. }
  350. }
  351. }
  352. template<typename eT>
  353. template<typename T1>
  354. inline
  355. void
  356. diagview<eT>::operator/=(const Base<eT,T1>& o)
  357. {
  358. arma_extra_debug_sigprint();
  359. diagview<eT>& d = *this;
  360. Mat<eT>& d_m = const_cast< Mat<eT>& >(d.m);
  361. const uword d_n_elem = d.n_elem;
  362. const uword d_row_offset = d.row_offset;
  363. const uword d_col_offset = d.col_offset;
  364. const Proxy<T1> P( o.get_ref() );
  365. arma_debug_check
  366. (
  367. ( (d_n_elem != P.get_n_elem()) || ((P.get_n_rows() != 1) && (P.get_n_cols() != 1)) ),
  368. "diagview: given object has incompatible size"
  369. );
  370. const bool is_alias = P.is_alias(d_m);
  371. if(is_alias) { arma_extra_debug_print("aliasing detected"); }
  372. if( (is_Mat<typename Proxy<T1>::stored_type>::value) || (Proxy<T1>::use_at) || (is_alias) )
  373. {
  374. const unwrap_check<typename Proxy<T1>::stored_type> tmp(P.Q, is_alias);
  375. const Mat<eT>& x = tmp.M;
  376. const eT* x_mem = x.memptr();
  377. uword ii,jj;
  378. for(ii=0, jj=1; jj < d_n_elem; ii+=2, jj+=2)
  379. {
  380. const eT tmp_i = x_mem[ii];
  381. const eT tmp_j = x_mem[jj];
  382. d_m.at( ii + d_row_offset, ii + d_col_offset) /= tmp_i;
  383. d_m.at( jj + d_row_offset, jj + d_col_offset) /= tmp_j;
  384. }
  385. if(ii < d_n_elem)
  386. {
  387. d_m.at( ii + d_row_offset, ii + d_col_offset) /= x_mem[ii];
  388. }
  389. }
  390. else
  391. {
  392. typename Proxy<T1>::ea_type Pea = P.get_ea();
  393. uword ii,jj;
  394. for(ii=0, jj=1; jj < d_n_elem; ii+=2, jj+=2)
  395. {
  396. const eT tmp_i = Pea[ii];
  397. const eT tmp_j = Pea[jj];
  398. d_m.at( ii + d_row_offset, ii + d_col_offset) /= tmp_i;
  399. d_m.at( jj + d_row_offset, jj + d_col_offset) /= tmp_j;
  400. }
  401. if(ii < d_n_elem)
  402. {
  403. d_m.at( ii + d_row_offset, ii + d_col_offset) /= Pea[ii];
  404. }
  405. }
  406. }
  407. //! extract a diagonal and store it as a column vector
  408. template<typename eT>
  409. inline
  410. void
  411. diagview<eT>::extract(Mat<eT>& out, const diagview<eT>& in)
  412. {
  413. arma_extra_debug_sigprint();
  414. // NOTE: we're assuming that the matrix has already been set to the correct size and there is no aliasing;
  415. // size setting and alias checking is done by either the Mat contructor or operator=()
  416. const Mat<eT>& in_m = in.m;
  417. const uword in_n_elem = in.n_elem;
  418. const uword in_row_offset = in.row_offset;
  419. const uword in_col_offset = in.col_offset;
  420. eT* out_mem = out.memptr();
  421. uword i,j;
  422. for(i=0, j=1; j < in_n_elem; i+=2, j+=2)
  423. {
  424. const eT tmp_i = in_m.at( i + in_row_offset, i + in_col_offset );
  425. const eT tmp_j = in_m.at( j + in_row_offset, j + in_col_offset );
  426. out_mem[i] = tmp_i;
  427. out_mem[j] = tmp_j;
  428. }
  429. if(i < in_n_elem)
  430. {
  431. out_mem[i] = in_m.at( i + in_row_offset, i + in_col_offset );
  432. }
  433. }
  434. //! X += Y.diag()
  435. template<typename eT>
  436. inline
  437. void
  438. diagview<eT>::plus_inplace(Mat<eT>& out, const diagview<eT>& in)
  439. {
  440. arma_extra_debug_sigprint();
  441. arma_debug_assert_same_size(out.n_rows, out.n_cols, in.n_rows, in.n_cols, "addition");
  442. const Mat<eT>& in_m = in.m;
  443. const uword in_n_elem = in.n_elem;
  444. const uword in_row_offset = in.row_offset;
  445. const uword in_col_offset = in.col_offset;
  446. eT* out_mem = out.memptr();
  447. uword i,j;
  448. for(i=0, j=1; j < in_n_elem; i+=2, j+=2)
  449. {
  450. const eT tmp_i = in_m.at( i + in_row_offset, i + in_col_offset );
  451. const eT tmp_j = in_m.at( j + in_row_offset, j + in_col_offset );
  452. out_mem[i] += tmp_i;
  453. out_mem[j] += tmp_j;
  454. }
  455. if(i < in_n_elem)
  456. {
  457. out_mem[i] += in_m.at( i + in_row_offset, i + in_col_offset );
  458. }
  459. }
  460. //! X -= Y.diag()
  461. template<typename eT>
  462. inline
  463. void
  464. diagview<eT>::minus_inplace(Mat<eT>& out, const diagview<eT>& in)
  465. {
  466. arma_extra_debug_sigprint();
  467. arma_debug_assert_same_size(out.n_rows, out.n_cols, in.n_rows, in.n_cols, "subtraction");
  468. const Mat<eT>& in_m = in.m;
  469. const uword in_n_elem = in.n_elem;
  470. const uword in_row_offset = in.row_offset;
  471. const uword in_col_offset = in.col_offset;
  472. eT* out_mem = out.memptr();
  473. uword i,j;
  474. for(i=0, j=1; j < in_n_elem; i+=2, j+=2)
  475. {
  476. const eT tmp_i = in_m.at( i + in_row_offset, i + in_col_offset );
  477. const eT tmp_j = in_m.at( j + in_row_offset, j + in_col_offset );
  478. out_mem[i] -= tmp_i;
  479. out_mem[j] -= tmp_j;
  480. }
  481. if(i < in_n_elem)
  482. {
  483. out_mem[i] -= in_m.at( i + in_row_offset, i + in_col_offset );
  484. }
  485. }
  486. //! X %= Y.diag()
  487. template<typename eT>
  488. inline
  489. void
  490. diagview<eT>::schur_inplace(Mat<eT>& out, const diagview<eT>& in)
  491. {
  492. arma_extra_debug_sigprint();
  493. arma_debug_assert_same_size(out.n_rows, out.n_cols, in.n_rows, in.n_cols, "element-wise multiplication");
  494. const Mat<eT>& in_m = in.m;
  495. const uword in_n_elem = in.n_elem;
  496. const uword in_row_offset = in.row_offset;
  497. const uword in_col_offset = in.col_offset;
  498. eT* out_mem = out.memptr();
  499. uword i,j;
  500. for(i=0, j=1; j < in_n_elem; i+=2, j+=2)
  501. {
  502. const eT tmp_i = in_m.at( i + in_row_offset, i + in_col_offset );
  503. const eT tmp_j = in_m.at( j + in_row_offset, j + in_col_offset );
  504. out_mem[i] *= tmp_i;
  505. out_mem[j] *= tmp_j;
  506. }
  507. if(i < in_n_elem)
  508. {
  509. out_mem[i] *= in_m.at( i + in_row_offset, i + in_col_offset );
  510. }
  511. }
  512. //! X /= Y.diag()
  513. template<typename eT>
  514. inline
  515. void
  516. diagview<eT>::div_inplace(Mat<eT>& out, const diagview<eT>& in)
  517. {
  518. arma_extra_debug_sigprint();
  519. arma_debug_assert_same_size(out.n_rows, out.n_cols, in.n_rows, in.n_cols, "element-wise division");
  520. const Mat<eT>& in_m = in.m;
  521. const uword in_n_elem = in.n_elem;
  522. const uword in_row_offset = in.row_offset;
  523. const uword in_col_offset = in.col_offset;
  524. eT* out_mem = out.memptr();
  525. uword i,j;
  526. for(i=0, j=1; j < in_n_elem; i+=2, j+=2)
  527. {
  528. const eT tmp_i = in_m.at( i + in_row_offset, i + in_col_offset );
  529. const eT tmp_j = in_m.at( j + in_row_offset, j + in_col_offset );
  530. out_mem[i] /= tmp_i;
  531. out_mem[j] /= tmp_j;
  532. }
  533. if(i < in_n_elem)
  534. {
  535. out_mem[i] /= in_m.at( i + in_row_offset, i + in_col_offset );
  536. }
  537. }
  538. template<typename eT>
  539. arma_inline
  540. eT
  541. diagview<eT>::at_alt(const uword ii) const
  542. {
  543. return m.at(ii+row_offset, ii+col_offset);
  544. }
  545. template<typename eT>
  546. arma_inline
  547. eT&
  548. diagview<eT>::operator[](const uword ii)
  549. {
  550. return (const_cast< Mat<eT>& >(m)).at(ii+row_offset, ii+col_offset);
  551. }
  552. template<typename eT>
  553. arma_inline
  554. eT
  555. diagview<eT>::operator[](const uword ii) const
  556. {
  557. return m.at(ii+row_offset, ii+col_offset);
  558. }
  559. template<typename eT>
  560. arma_inline
  561. eT&
  562. diagview<eT>::at(const uword ii)
  563. {
  564. return (const_cast< Mat<eT>& >(m)).at(ii+row_offset, ii+col_offset);
  565. }
  566. template<typename eT>
  567. arma_inline
  568. eT
  569. diagview<eT>::at(const uword ii) const
  570. {
  571. return m.at(ii+row_offset, ii+col_offset);
  572. }
  573. template<typename eT>
  574. arma_inline
  575. eT&
  576. diagview<eT>::operator()(const uword ii)
  577. {
  578. arma_debug_check( (ii >= n_elem), "diagview::operator(): out of bounds" );
  579. return (const_cast< Mat<eT>& >(m)).at(ii+row_offset, ii+col_offset);
  580. }
  581. template<typename eT>
  582. arma_inline
  583. eT
  584. diagview<eT>::operator()(const uword ii) const
  585. {
  586. arma_debug_check( (ii >= n_elem), "diagview::operator(): out of bounds" );
  587. return m.at(ii+row_offset, ii+col_offset);
  588. }
  589. template<typename eT>
  590. arma_inline
  591. eT&
  592. diagview<eT>::at(const uword row, const uword)
  593. {
  594. return (const_cast< Mat<eT>& >(m)).at(row+row_offset, row+col_offset);
  595. }
  596. template<typename eT>
  597. arma_inline
  598. eT
  599. diagview<eT>::at(const uword row, const uword) const
  600. {
  601. return m.at(row+row_offset, row+col_offset);
  602. }
  603. template<typename eT>
  604. arma_inline
  605. eT&
  606. diagview<eT>::operator()(const uword row, const uword col)
  607. {
  608. arma_debug_check( ((row >= n_elem) || (col > 0)), "diagview::operator(): out of bounds" );
  609. return (const_cast< Mat<eT>& >(m)).at(row+row_offset, row+col_offset);
  610. }
  611. template<typename eT>
  612. arma_inline
  613. eT
  614. diagview<eT>::operator()(const uword row, const uword col) const
  615. {
  616. arma_debug_check( ((row >= n_elem) || (col > 0)), "diagview::operator(): out of bounds" );
  617. return m.at(row+row_offset, row+col_offset);
  618. }
  619. template<typename eT>
  620. arma_inline
  621. const Op<diagview<eT>,op_htrans>
  622. diagview<eT>::t() const
  623. {
  624. return Op<diagview<eT>,op_htrans>(*this);
  625. }
  626. template<typename eT>
  627. arma_inline
  628. const Op<diagview<eT>,op_htrans>
  629. diagview<eT>::ht() const
  630. {
  631. return Op<diagview<eT>,op_htrans>(*this);
  632. }
  633. template<typename eT>
  634. arma_inline
  635. const Op<diagview<eT>,op_strans>
  636. diagview<eT>::st() const
  637. {
  638. return Op<diagview<eT>,op_strans>(*this);
  639. }
  640. template<typename eT>
  641. inline
  642. void
  643. diagview<eT>::replace(const eT old_val, const eT new_val)
  644. {
  645. arma_extra_debug_sigprint();
  646. Mat<eT>& x = const_cast< Mat<eT>& >(m);
  647. const uword local_n_elem = n_elem;
  648. if(arma_isnan(old_val))
  649. {
  650. for(uword ii=0; ii < local_n_elem; ++ii)
  651. {
  652. eT& val = x.at(ii+row_offset, ii+col_offset);
  653. val = (arma_isnan(val)) ? new_val : val;
  654. }
  655. }
  656. else
  657. {
  658. for(uword ii=0; ii < local_n_elem; ++ii)
  659. {
  660. eT& val = x.at(ii+row_offset, ii+col_offset);
  661. val = (val == old_val) ? new_val : val;
  662. }
  663. }
  664. }
  665. template<typename eT>
  666. inline
  667. void
  668. diagview<eT>::fill(const eT val)
  669. {
  670. arma_extra_debug_sigprint();
  671. Mat<eT>& x = const_cast< Mat<eT>& >(m);
  672. const uword local_n_elem = n_elem;
  673. for(uword ii=0; ii < local_n_elem; ++ii)
  674. {
  675. x.at(ii+row_offset, ii+col_offset) = val;
  676. }
  677. }
  678. template<typename eT>
  679. inline
  680. void
  681. diagview<eT>::zeros()
  682. {
  683. arma_extra_debug_sigprint();
  684. (*this).fill(eT(0));
  685. }
  686. template<typename eT>
  687. inline
  688. void
  689. diagview<eT>::ones()
  690. {
  691. arma_extra_debug_sigprint();
  692. (*this).fill(eT(1));
  693. }
  694. template<typename eT>
  695. inline
  696. void
  697. diagview<eT>::randu()
  698. {
  699. arma_extra_debug_sigprint();
  700. Mat<eT>& x = const_cast< Mat<eT>& >(m);
  701. const uword local_n_elem = n_elem;
  702. for(uword ii=0; ii < local_n_elem; ++ii)
  703. {
  704. x.at(ii+row_offset, ii+col_offset) = eT(arma_rng::randu<eT>());
  705. }
  706. }
  707. template<typename eT>
  708. inline
  709. void
  710. diagview<eT>::randn()
  711. {
  712. arma_extra_debug_sigprint();
  713. Mat<eT>& x = const_cast< Mat<eT>& >(m);
  714. const uword local_n_elem = n_elem;
  715. for(uword ii=0; ii < local_n_elem; ++ii)
  716. {
  717. x.at(ii+row_offset, ii+col_offset) = eT(arma_rng::randn<eT>());
  718. }
  719. }
  720. //! @}