1 // SPDX-License-Identifier: Apache-2.0
2 //
3 // Copyright 2008-2016 Conrad Sanderson (http://conradsanderson.id.au)
4 // Copyright 2008-2016 National ICT Australia (NICTA)
5 //
6 // Licensed under the Apache License, Version 2.0 (the "License");
7 // you may not use this file except in compliance with the License.
8 // You may obtain a copy of the License at
9 // http://www.apache.org/licenses/LICENSE-2.0
10 //
11 // Unless required by applicable law or agreed to in writing, software
12 // distributed under the License is distributed on an "AS IS" BASIS,
13 // WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14 // See the License for the specific language governing permissions and
15 // limitations under the License.
16 // ------------------------------------------------------------------------
17 
18 
19 //! \addtogroup diagview
20 //! @{
21 
22 
23 template<typename eT>
24 inline
~diagview()25 diagview<eT>::~diagview()
26   {
27   arma_extra_debug_sigprint_this(this);
28   }
29 
30 
31 
32 template<typename eT>
33 arma_inline
diagview(const Mat<eT> & in_m,const uword in_row_offset,const uword in_col_offset,const uword in_len)34 diagview<eT>::diagview(const Mat<eT>& in_m, const uword in_row_offset, const uword in_col_offset, const uword in_len)
35   : m         (in_m         )
36   , row_offset(in_row_offset)
37   , col_offset(in_col_offset)
38   , n_rows    (in_len       )
39   , n_elem    (in_len       )
40   {
41   arma_extra_debug_sigprint_this(this);
42   }
43 
44 
45 
46 template<typename eT>
47 inline
diagview(const diagview<eT> & in)48 diagview<eT>::diagview(const diagview<eT>& in)
49   : m         (in.m         )
50   , row_offset(in.row_offset)
51   , col_offset(in.col_offset)
52   , n_rows    (in.n_rows    )
53   , n_elem    (in.n_elem    )
54   {
55   arma_extra_debug_sigprint(arma_str::format("this = %x   in = %x") % this % &in);
56   }
57 
58 
59 
60 template<typename eT>
61 inline
diagview(diagview<eT> && in)62 diagview<eT>::diagview(diagview<eT>&& in)
63   : m         (in.m         )
64   , row_offset(in.row_offset)
65   , col_offset(in.col_offset)
66   , n_rows    (in.n_rows    )
67   , n_elem    (in.n_elem    )
68   {
69   arma_extra_debug_sigprint(arma_str::format("this = %x   in = %x") % this % &in);
70 
71   // for paranoia
72 
73   access::rw(in.row_offset) = 0;
74   access::rw(in.col_offset) = 0;
75   access::rw(in.n_rows    ) = 0;
76   access::rw(in.n_elem    ) = 0;
77   }
78 
79 
80 
81 //! set a diagonal of our matrix using a diagonal from a foreign matrix
82 template<typename eT>
83 inline
84 void
operator =(const diagview<eT> & x)85 diagview<eT>::operator= (const diagview<eT>& x)
86   {
87   arma_extra_debug_sigprint();
88 
89   diagview<eT>& d = *this;
90 
91   arma_debug_check( (d.n_elem != x.n_elem), "diagview: diagonals have incompatible lengths" );
92 
93         Mat<eT>& d_m = const_cast< Mat<eT>& >(d.m);
94   const Mat<eT>& x_m = x.m;
95 
96   if(&d_m != &x_m)
97     {
98     const uword d_n_elem     = d.n_elem;
99     const uword d_row_offset = d.row_offset;
100     const uword d_col_offset = d.col_offset;
101 
102     const uword x_row_offset = x.row_offset;
103     const uword x_col_offset = x.col_offset;
104 
105     uword ii,jj;
106     for(ii=0, jj=1; jj < d_n_elem; ii+=2, jj+=2)
107       {
108       const eT tmp_i = x_m.at(ii + x_row_offset, ii + x_col_offset);
109       const eT tmp_j = x_m.at(jj + x_row_offset, jj + x_col_offset);
110 
111       d_m.at(ii + d_row_offset, ii + d_col_offset) = tmp_i;
112       d_m.at(jj + d_row_offset, jj + d_col_offset) = tmp_j;
113       }
114 
115     if(ii < d_n_elem)
116       {
117       d_m.at(ii + d_row_offset, ii + d_col_offset) = x_m.at(ii + x_row_offset, ii + x_col_offset);
118       }
119     }
120   else
121     {
122     const Mat<eT> tmp = x;
123 
124     (*this).operator=(tmp);
125     }
126   }
127 
128 
129 
130 template<typename eT>
131 inline
132 void
operator +=(const eT val)133 diagview<eT>::operator+=(const eT val)
134   {
135   arma_extra_debug_sigprint();
136 
137   Mat<eT>& t_m = const_cast< Mat<eT>& >(m);
138 
139   const uword t_n_elem     = n_elem;
140   const uword t_row_offset = row_offset;
141   const uword t_col_offset = col_offset;
142 
143   for(uword ii=0; ii < t_n_elem; ++ii)
144     {
145     t_m.at( ii + t_row_offset,  ii + t_col_offset) += val;
146     }
147   }
148 
149 
150 
151 template<typename eT>
152 inline
153 void
operator -=(const eT val)154 diagview<eT>::operator-=(const eT val)
155   {
156   arma_extra_debug_sigprint();
157 
158   Mat<eT>& t_m = const_cast< Mat<eT>& >(m);
159 
160   const uword t_n_elem     = n_elem;
161   const uword t_row_offset = row_offset;
162   const uword t_col_offset = col_offset;
163 
164   for(uword ii=0; ii < t_n_elem; ++ii)
165     {
166     t_m.at( ii + t_row_offset,  ii + t_col_offset) -= val;
167     }
168   }
169 
170 
171 
172 template<typename eT>
173 inline
174 void
operator *=(const eT val)175 diagview<eT>::operator*=(const eT val)
176   {
177   arma_extra_debug_sigprint();
178 
179   Mat<eT>& t_m = const_cast< Mat<eT>& >(m);
180 
181   const uword t_n_elem     = n_elem;
182   const uword t_row_offset = row_offset;
183   const uword t_col_offset = col_offset;
184 
185   for(uword ii=0; ii < t_n_elem; ++ii)
186     {
187     t_m.at( ii + t_row_offset,  ii + t_col_offset) *= val;
188     }
189   }
190 
191 
192 
193 template<typename eT>
194 inline
195 void
operator /=(const eT val)196 diagview<eT>::operator/=(const eT val)
197   {
198   arma_extra_debug_sigprint();
199 
200   Mat<eT>& t_m = const_cast< Mat<eT>& >(m);
201 
202   const uword t_n_elem     = n_elem;
203   const uword t_row_offset = row_offset;
204   const uword t_col_offset = col_offset;
205 
206   for(uword ii=0; ii < t_n_elem; ++ii)
207     {
208     t_m.at( ii + t_row_offset,  ii + t_col_offset) /= val;
209     }
210   }
211 
212 
213 
214 //! set a diagonal of our matrix using data from a foreign object
215 template<typename eT>
216 template<typename T1>
217 inline
218 void
operator =(const Base<eT,T1> & o)219 diagview<eT>::operator= (const Base<eT,T1>& o)
220   {
221   arma_extra_debug_sigprint();
222 
223   diagview<eT>& d = *this;
224 
225   Mat<eT>& d_m = const_cast< Mat<eT>& >(d.m);
226 
227   const uword d_n_elem     = d.n_elem;
228   const uword d_row_offset = d.row_offset;
229   const uword d_col_offset = d.col_offset;
230 
231   const Proxy<T1> P( o.get_ref() );
232 
233   arma_debug_check
234     (
235     ( (d_n_elem != P.get_n_elem()) || ((P.get_n_rows() != 1) && (P.get_n_cols() != 1)) ),
236     "diagview: given object has incompatible size"
237     );
238 
239   const bool is_alias = P.is_alias(d_m);
240 
241   if(is_alias)  { arma_extra_debug_print("aliasing detected"); }
242 
243   if( (is_Mat<typename Proxy<T1>::stored_type>::value) || (Proxy<T1>::use_at) || (is_alias) )
244     {
245     const unwrap_check<typename Proxy<T1>::stored_type> tmp(P.Q, is_alias);
246     const Mat<eT>& x = tmp.M;
247 
248     const eT* x_mem = x.memptr();
249 
250     uword ii,jj;
251     for(ii=0, jj=1; jj < d_n_elem; ii+=2, jj+=2)
252       {
253       const eT tmp_i = x_mem[ii];
254       const eT tmp_j = x_mem[jj];
255 
256       d_m.at( ii + d_row_offset,  ii + d_col_offset) = tmp_i;
257       d_m.at( jj + d_row_offset,  jj + d_col_offset) = tmp_j;
258       }
259 
260     if(ii < d_n_elem)
261       {
262       d_m.at( ii + d_row_offset,  ii + d_col_offset) = x_mem[ii];
263       }
264     }
265   else
266     {
267     typename Proxy<T1>::ea_type Pea = P.get_ea();
268 
269     uword ii,jj;
270     for(ii=0, jj=1; jj < d_n_elem; ii+=2, jj+=2)
271       {
272       const eT tmp_i = Pea[ii];
273       const eT tmp_j = Pea[jj];
274 
275       d_m.at( ii + d_row_offset,  ii + d_col_offset) = tmp_i;
276       d_m.at( jj + d_row_offset,  jj + d_col_offset) = tmp_j;
277       }
278 
279     if(ii < d_n_elem)
280       {
281       d_m.at( ii + d_row_offset,  ii + d_col_offset) = Pea[ii];
282       }
283     }
284   }
285 
286 
287 
288 template<typename eT>
289 template<typename T1>
290 inline
291 void
operator +=(const Base<eT,T1> & o)292 diagview<eT>::operator+=(const Base<eT,T1>& o)
293   {
294   arma_extra_debug_sigprint();
295 
296   diagview<eT>& d = *this;
297 
298   Mat<eT>& d_m = const_cast< Mat<eT>& >(d.m);
299 
300   const uword d_n_elem     = d.n_elem;
301   const uword d_row_offset = d.row_offset;
302   const uword d_col_offset = d.col_offset;
303 
304   const Proxy<T1> P( o.get_ref() );
305 
306   arma_debug_check
307     (
308     ( (d_n_elem != P.get_n_elem()) || ((P.get_n_rows() != 1) && (P.get_n_cols() != 1)) ),
309     "diagview: given object has incompatible size"
310     );
311 
312   const bool is_alias = P.is_alias(d_m);
313 
314   if(is_alias)  { arma_extra_debug_print("aliasing detected"); }
315 
316   if( (is_Mat<typename Proxy<T1>::stored_type>::value) || (Proxy<T1>::use_at) || (is_alias) )
317     {
318     const unwrap_check<typename Proxy<T1>::stored_type> tmp(P.Q, is_alias);
319     const Mat<eT>& x = tmp.M;
320 
321     const eT* x_mem = x.memptr();
322 
323     uword ii,jj;
324     for(ii=0, jj=1; jj < d_n_elem; ii+=2, jj+=2)
325       {
326       const eT tmp_i = x_mem[ii];
327       const eT tmp_j = x_mem[jj];
328 
329       d_m.at( ii + d_row_offset,  ii + d_col_offset) += tmp_i;
330       d_m.at( jj + d_row_offset,  jj + d_col_offset) += tmp_j;
331       }
332 
333     if(ii < d_n_elem)
334       {
335       d_m.at( ii + d_row_offset,  ii + d_col_offset) += x_mem[ii];
336       }
337     }
338   else
339     {
340     typename Proxy<T1>::ea_type Pea = P.get_ea();
341 
342     uword ii,jj;
343     for(ii=0, jj=1; jj < d_n_elem; ii+=2, jj+=2)
344       {
345       const eT tmp_i = Pea[ii];
346       const eT tmp_j = Pea[jj];
347 
348       d_m.at( ii + d_row_offset,  ii + d_col_offset) += tmp_i;
349       d_m.at( jj + d_row_offset,  jj + d_col_offset) += tmp_j;
350       }
351 
352     if(ii < d_n_elem)
353       {
354       d_m.at( ii + d_row_offset,  ii + d_col_offset) += Pea[ii];
355       }
356     }
357   }
358 
359 
360 
361 template<typename eT>
362 template<typename T1>
363 inline
364 void
operator -=(const Base<eT,T1> & o)365 diagview<eT>::operator-=(const Base<eT,T1>& o)
366   {
367   arma_extra_debug_sigprint();
368 
369   diagview<eT>& d = *this;
370 
371   Mat<eT>& d_m = const_cast< Mat<eT>& >(d.m);
372 
373   const uword d_n_elem     = d.n_elem;
374   const uword d_row_offset = d.row_offset;
375   const uword d_col_offset = d.col_offset;
376 
377   const Proxy<T1> P( o.get_ref() );
378 
379   arma_debug_check
380     (
381     ( (d_n_elem != P.get_n_elem()) || ((P.get_n_rows() != 1) && (P.get_n_cols() != 1)) ),
382     "diagview: given object has incompatible size"
383     );
384 
385   const bool is_alias = P.is_alias(d_m);
386 
387   if(is_alias)  { arma_extra_debug_print("aliasing detected"); }
388 
389   if( (is_Mat<typename Proxy<T1>::stored_type>::value) || (Proxy<T1>::use_at) || (is_alias) )
390     {
391     const unwrap_check<typename Proxy<T1>::stored_type> tmp(P.Q, is_alias);
392     const Mat<eT>& x = tmp.M;
393 
394     const eT* x_mem = x.memptr();
395 
396     uword ii,jj;
397     for(ii=0, jj=1; jj < d_n_elem; ii+=2, jj+=2)
398       {
399       const eT tmp_i = x_mem[ii];
400       const eT tmp_j = x_mem[jj];
401 
402       d_m.at( ii + d_row_offset,  ii + d_col_offset) -= tmp_i;
403       d_m.at( jj + d_row_offset,  jj + d_col_offset) -= tmp_j;
404       }
405 
406     if(ii < d_n_elem)
407       {
408       d_m.at( ii + d_row_offset,  ii + d_col_offset) -= x_mem[ii];
409       }
410     }
411   else
412     {
413     typename Proxy<T1>::ea_type Pea = P.get_ea();
414 
415     uword ii,jj;
416     for(ii=0, jj=1; jj < d_n_elem; ii+=2, jj+=2)
417       {
418       const eT tmp_i = Pea[ii];
419       const eT tmp_j = Pea[jj];
420 
421       d_m.at( ii + d_row_offset,  ii + d_col_offset) -= tmp_i;
422       d_m.at( jj + d_row_offset,  jj + d_col_offset) -= tmp_j;
423       }
424 
425     if(ii < d_n_elem)
426       {
427       d_m.at( ii + d_row_offset,  ii + d_col_offset) -= Pea[ii];
428       }
429     }
430   }
431 
432 
433 
434 template<typename eT>
435 template<typename T1>
436 inline
437 void
operator %=(const Base<eT,T1> & o)438 diagview<eT>::operator%=(const Base<eT,T1>& o)
439   {
440   arma_extra_debug_sigprint();
441 
442   diagview<eT>& d = *this;
443 
444   Mat<eT>& d_m = const_cast< Mat<eT>& >(d.m);
445 
446   const uword d_n_elem     = d.n_elem;
447   const uword d_row_offset = d.row_offset;
448   const uword d_col_offset = d.col_offset;
449 
450   const Proxy<T1> P( o.get_ref() );
451 
452   arma_debug_check
453     (
454     ( (d_n_elem != P.get_n_elem()) || ((P.get_n_rows() != 1) && (P.get_n_cols() != 1)) ),
455     "diagview: given object has incompatible size"
456     );
457 
458   const bool is_alias = P.is_alias(d_m);
459 
460   if(is_alias)  { arma_extra_debug_print("aliasing detected"); }
461 
462   if( (is_Mat<typename Proxy<T1>::stored_type>::value) || (Proxy<T1>::use_at) || (is_alias) )
463     {
464     const unwrap_check<typename Proxy<T1>::stored_type> tmp(P.Q, is_alias);
465     const Mat<eT>& x = tmp.M;
466 
467     const eT* x_mem = x.memptr();
468 
469     uword ii,jj;
470     for(ii=0, jj=1; jj < d_n_elem; ii+=2, jj+=2)
471       {
472       const eT tmp_i = x_mem[ii];
473       const eT tmp_j = x_mem[jj];
474 
475       d_m.at( ii + d_row_offset,  ii + d_col_offset) *= tmp_i;
476       d_m.at( jj + d_row_offset,  jj + d_col_offset) *= tmp_j;
477       }
478 
479     if(ii < d_n_elem)
480       {
481       d_m.at( ii + d_row_offset,  ii + d_col_offset) *= x_mem[ii];
482       }
483     }
484   else
485     {
486     typename Proxy<T1>::ea_type Pea = P.get_ea();
487 
488     uword ii,jj;
489     for(ii=0, jj=1; jj < d_n_elem; ii+=2, jj+=2)
490       {
491       const eT tmp_i = Pea[ii];
492       const eT tmp_j = Pea[jj];
493 
494       d_m.at( ii + d_row_offset,  ii + d_col_offset) *= tmp_i;
495       d_m.at( jj + d_row_offset,  jj + d_col_offset) *= tmp_j;
496       }
497 
498     if(ii < d_n_elem)
499       {
500       d_m.at( ii + d_row_offset,  ii + d_col_offset) *= Pea[ii];
501       }
502     }
503   }
504 
505 
506 
507 template<typename eT>
508 template<typename T1>
509 inline
510 void
operator /=(const Base<eT,T1> & o)511 diagview<eT>::operator/=(const Base<eT,T1>& o)
512   {
513   arma_extra_debug_sigprint();
514 
515   diagview<eT>& d = *this;
516 
517   Mat<eT>& d_m = const_cast< Mat<eT>& >(d.m);
518 
519   const uword d_n_elem     = d.n_elem;
520   const uword d_row_offset = d.row_offset;
521   const uword d_col_offset = d.col_offset;
522 
523   const Proxy<T1> P( o.get_ref() );
524 
525   arma_debug_check
526     (
527     ( (d_n_elem != P.get_n_elem()) || ((P.get_n_rows() != 1) && (P.get_n_cols() != 1)) ),
528     "diagview: given object has incompatible size"
529     );
530 
531   const bool is_alias = P.is_alias(d_m);
532 
533   if(is_alias)  { arma_extra_debug_print("aliasing detected"); }
534 
535   if( (is_Mat<typename Proxy<T1>::stored_type>::value) || (Proxy<T1>::use_at) || (is_alias) )
536     {
537     const unwrap_check<typename Proxy<T1>::stored_type> tmp(P.Q, is_alias);
538     const Mat<eT>& x = tmp.M;
539 
540     const eT* x_mem = x.memptr();
541 
542     uword ii,jj;
543     for(ii=0, jj=1; jj < d_n_elem; ii+=2, jj+=2)
544       {
545       const eT tmp_i = x_mem[ii];
546       const eT tmp_j = x_mem[jj];
547 
548       d_m.at( ii + d_row_offset,  ii + d_col_offset) /= tmp_i;
549       d_m.at( jj + d_row_offset,  jj + d_col_offset) /= tmp_j;
550       }
551 
552     if(ii < d_n_elem)
553       {
554       d_m.at( ii + d_row_offset,  ii + d_col_offset) /= x_mem[ii];
555       }
556     }
557   else
558     {
559     typename Proxy<T1>::ea_type Pea = P.get_ea();
560 
561     uword ii,jj;
562     for(ii=0, jj=1; jj < d_n_elem; ii+=2, jj+=2)
563       {
564       const eT tmp_i = Pea[ii];
565       const eT tmp_j = Pea[jj];
566 
567       d_m.at( ii + d_row_offset,  ii + d_col_offset) /= tmp_i;
568       d_m.at( jj + d_row_offset,  jj + d_col_offset) /= tmp_j;
569       }
570 
571     if(ii < d_n_elem)
572       {
573       d_m.at( ii + d_row_offset,  ii + d_col_offset) /= Pea[ii];
574       }
575     }
576   }
577 
578 
579 
580 //! extract a diagonal and store it as a column vector
581 template<typename eT>
582 inline
583 void
extract(Mat<eT> & out,const diagview<eT> & in)584 diagview<eT>::extract(Mat<eT>& out, const diagview<eT>& in)
585   {
586   arma_extra_debug_sigprint();
587 
588   // NOTE: we're assuming that the matrix has already been set to the correct size and there is no aliasing;
589   // size setting and alias checking is done by either the Mat contructor or operator=()
590 
591   const Mat<eT>& in_m = in.m;
592 
593   const uword in_n_elem     = in.n_elem;
594   const uword in_row_offset = in.row_offset;
595   const uword in_col_offset = in.col_offset;
596 
597   eT* out_mem = out.memptr();
598 
599   uword i,j;
600   for(i=0, j=1; j < in_n_elem; i+=2, j+=2)
601     {
602     const eT tmp_i = in_m.at( i + in_row_offset, i + in_col_offset );
603     const eT tmp_j = in_m.at( j + in_row_offset, j + in_col_offset );
604 
605     out_mem[i] = tmp_i;
606     out_mem[j] = tmp_j;
607     }
608 
609   if(i < in_n_elem)
610     {
611     out_mem[i] = in_m.at( i + in_row_offset, i + in_col_offset );
612     }
613   }
614 
615 
616 
617 //! X += Y.diag()
618 template<typename eT>
619 inline
620 void
plus_inplace(Mat<eT> & out,const diagview<eT> & in)621 diagview<eT>::plus_inplace(Mat<eT>& out, const diagview<eT>& in)
622   {
623   arma_extra_debug_sigprint();
624 
625   arma_debug_assert_same_size(out.n_rows, out.n_cols, in.n_rows, in.n_cols, "addition");
626 
627   const Mat<eT>& in_m = in.m;
628 
629   const uword in_n_elem     = in.n_elem;
630   const uword in_row_offset = in.row_offset;
631   const uword in_col_offset = in.col_offset;
632 
633   eT* out_mem = out.memptr();
634 
635   uword i,j;
636   for(i=0, j=1; j < in_n_elem; i+=2, j+=2)
637     {
638     const eT tmp_i = in_m.at( i + in_row_offset, i + in_col_offset );
639     const eT tmp_j = in_m.at( j + in_row_offset, j + in_col_offset );
640 
641     out_mem[i] += tmp_i;
642     out_mem[j] += tmp_j;
643     }
644 
645   if(i < in_n_elem)
646     {
647     out_mem[i] += in_m.at( i + in_row_offset, i + in_col_offset );
648     }
649   }
650 
651 
652 
653 //! X -= Y.diag()
654 template<typename eT>
655 inline
656 void
minus_inplace(Mat<eT> & out,const diagview<eT> & in)657 diagview<eT>::minus_inplace(Mat<eT>& out, const diagview<eT>& in)
658   {
659   arma_extra_debug_sigprint();
660 
661   arma_debug_assert_same_size(out.n_rows, out.n_cols, in.n_rows, in.n_cols, "subtraction");
662 
663   const Mat<eT>& in_m = in.m;
664 
665   const uword in_n_elem     = in.n_elem;
666   const uword in_row_offset = in.row_offset;
667   const uword in_col_offset = in.col_offset;
668 
669   eT* out_mem = out.memptr();
670 
671   uword i,j;
672   for(i=0, j=1; j < in_n_elem; i+=2, j+=2)
673     {
674     const eT tmp_i = in_m.at( i + in_row_offset, i + in_col_offset );
675     const eT tmp_j = in_m.at( j + in_row_offset, j + in_col_offset );
676 
677     out_mem[i] -= tmp_i;
678     out_mem[j] -= tmp_j;
679     }
680 
681   if(i < in_n_elem)
682     {
683     out_mem[i] -= in_m.at( i + in_row_offset, i + in_col_offset );
684     }
685   }
686 
687 
688 
689 //! X %= Y.diag()
690 template<typename eT>
691 inline
692 void
schur_inplace(Mat<eT> & out,const diagview<eT> & in)693 diagview<eT>::schur_inplace(Mat<eT>& out, const diagview<eT>& in)
694   {
695   arma_extra_debug_sigprint();
696 
697   arma_debug_assert_same_size(out.n_rows, out.n_cols, in.n_rows, in.n_cols, "element-wise multiplication");
698 
699   const Mat<eT>& in_m = in.m;
700 
701   const uword in_n_elem     = in.n_elem;
702   const uword in_row_offset = in.row_offset;
703   const uword in_col_offset = in.col_offset;
704 
705   eT* out_mem = out.memptr();
706 
707   uword i,j;
708   for(i=0, j=1; j < in_n_elem; i+=2, j+=2)
709     {
710     const eT tmp_i = in_m.at( i + in_row_offset, i + in_col_offset );
711     const eT tmp_j = in_m.at( j + in_row_offset, j + in_col_offset );
712 
713     out_mem[i] *= tmp_i;
714     out_mem[j] *= tmp_j;
715     }
716 
717   if(i < in_n_elem)
718     {
719     out_mem[i] *= in_m.at( i + in_row_offset, i + in_col_offset );
720     }
721   }
722 
723 
724 
725 //! X /= Y.diag()
726 template<typename eT>
727 inline
728 void
div_inplace(Mat<eT> & out,const diagview<eT> & in)729 diagview<eT>::div_inplace(Mat<eT>& out, const diagview<eT>& in)
730   {
731   arma_extra_debug_sigprint();
732 
733   arma_debug_assert_same_size(out.n_rows, out.n_cols, in.n_rows, in.n_cols, "element-wise division");
734 
735   const Mat<eT>& in_m = in.m;
736 
737   const uword in_n_elem     = in.n_elem;
738   const uword in_row_offset = in.row_offset;
739   const uword in_col_offset = in.col_offset;
740 
741   eT* out_mem = out.memptr();
742 
743   uword i,j;
744   for(i=0, j=1; j < in_n_elem; i+=2, j+=2)
745     {
746     const eT tmp_i = in_m.at( i + in_row_offset, i + in_col_offset );
747     const eT tmp_j = in_m.at( j + in_row_offset, j + in_col_offset );
748 
749     out_mem[i] /= tmp_i;
750     out_mem[j] /= tmp_j;
751     }
752 
753   if(i < in_n_elem)
754     {
755     out_mem[i] /= in_m.at( i + in_row_offset, i + in_col_offset );
756     }
757   }
758 
759 
760 
761 template<typename eT>
762 arma_inline
763 eT
at_alt(const uword ii) const764 diagview<eT>::at_alt(const uword ii) const
765   {
766   return m.at(ii+row_offset, ii+col_offset);
767   }
768 
769 
770 
771 template<typename eT>
772 arma_inline
773 eT&
operator [](const uword ii)774 diagview<eT>::operator[](const uword ii)
775   {
776   return (const_cast< Mat<eT>& >(m)).at(ii+row_offset, ii+col_offset);
777   }
778 
779 
780 
781 template<typename eT>
782 arma_inline
783 eT
operator [](const uword ii) const784 diagview<eT>::operator[](const uword ii) const
785   {
786   return m.at(ii+row_offset, ii+col_offset);
787   }
788 
789 
790 
791 template<typename eT>
792 arma_inline
793 eT&
at(const uword ii)794 diagview<eT>::at(const uword ii)
795   {
796   return (const_cast< Mat<eT>& >(m)).at(ii+row_offset, ii+col_offset);
797   }
798 
799 
800 
801 template<typename eT>
802 arma_inline
803 eT
at(const uword ii) const804 diagview<eT>::at(const uword ii) const
805   {
806   return m.at(ii+row_offset, ii+col_offset);
807   }
808 
809 
810 
811 template<typename eT>
812 arma_inline
813 eT&
operator ()(const uword ii)814 diagview<eT>::operator()(const uword ii)
815   {
816   arma_debug_check_bounds( (ii >= n_elem), "diagview::operator(): out of bounds" );
817 
818   return (const_cast< Mat<eT>& >(m)).at(ii+row_offset, ii+col_offset);
819   }
820 
821 
822 
823 template<typename eT>
824 arma_inline
825 eT
operator ()(const uword ii) const826 diagview<eT>::operator()(const uword ii) const
827   {
828   arma_debug_check_bounds( (ii >= n_elem), "diagview::operator(): out of bounds" );
829 
830   return m.at(ii+row_offset, ii+col_offset);
831   }
832 
833 
834 
835 template<typename eT>
836 arma_inline
837 eT&
at(const uword row,const uword)838 diagview<eT>::at(const uword row, const uword)
839   {
840   return (const_cast< Mat<eT>& >(m)).at(row+row_offset, row+col_offset);
841   }
842 
843 
844 
845 template<typename eT>
846 arma_inline
847 eT
at(const uword row,const uword) const848 diagview<eT>::at(const uword row, const uword) const
849   {
850   return m.at(row+row_offset, row+col_offset);
851   }
852 
853 
854 
855 template<typename eT>
856 arma_inline
857 eT&
operator ()(const uword row,const uword col)858 diagview<eT>::operator()(const uword row, const uword col)
859   {
860   arma_debug_check_bounds( ((row >= n_elem) || (col > 0)), "diagview::operator(): out of bounds" );
861 
862   return (const_cast< Mat<eT>& >(m)).at(row+row_offset, row+col_offset);
863   }
864 
865 
866 
867 template<typename eT>
868 arma_inline
869 eT
operator ()(const uword row,const uword col) const870 diagview<eT>::operator()(const uword row, const uword col) const
871   {
872   arma_debug_check_bounds( ((row >= n_elem) || (col > 0)), "diagview::operator(): out of bounds" );
873 
874   return m.at(row+row_offset, row+col_offset);
875   }
876 
877 
878 
879 template<typename eT>
880 inline
881 void
replace(const eT old_val,const eT new_val)882 diagview<eT>::replace(const eT old_val, const eT new_val)
883   {
884   arma_extra_debug_sigprint();
885 
886   Mat<eT>& x = const_cast< Mat<eT>& >(m);
887 
888   const uword local_n_elem = n_elem;
889 
890   if(arma_isnan(old_val))
891     {
892     for(uword ii=0; ii < local_n_elem; ++ii)
893       {
894       eT& val = x.at(ii+row_offset, ii+col_offset);
895 
896       val = (arma_isnan(val)) ? new_val : val;
897       }
898     }
899   else
900     {
901     for(uword ii=0; ii < local_n_elem; ++ii)
902       {
903       eT& val = x.at(ii+row_offset, ii+col_offset);
904 
905       val = (val == old_val) ? new_val : val;
906       }
907     }
908   }
909 
910 
911 
912 template<typename eT>
913 inline
914 void
clean(const typename get_pod_type<eT>::result threshold)915 diagview<eT>::clean(const typename get_pod_type<eT>::result threshold)
916   {
917   arma_extra_debug_sigprint();
918 
919   Mat<eT> tmp(*this);
920 
921   tmp.clean(threshold);
922 
923   (*this).operator=(tmp);
924   }
925 
926 
927 
928 template<typename eT>
929 inline
930 void
clamp(const eT min_val,const eT max_val)931 diagview<eT>::clamp(const eT min_val, const eT max_val)
932   {
933   arma_extra_debug_sigprint();
934 
935   Mat<eT> tmp(*this);
936 
937   tmp.clamp(min_val, max_val);
938 
939   (*this).operator=(tmp);
940   }
941 
942 
943 
944 template<typename eT>
945 inline
946 void
fill(const eT val)947 diagview<eT>::fill(const eT val)
948   {
949   arma_extra_debug_sigprint();
950 
951   Mat<eT>& x = const_cast< Mat<eT>& >(m);
952 
953   const uword local_n_elem = n_elem;
954 
955   for(uword ii=0; ii < local_n_elem; ++ii)
956     {
957     x.at(ii+row_offset, ii+col_offset) = val;
958     }
959   }
960 
961 
962 
963 template<typename eT>
964 inline
965 void
zeros()966 diagview<eT>::zeros()
967   {
968   arma_extra_debug_sigprint();
969 
970   (*this).fill(eT(0));
971   }
972 
973 
974 
975 template<typename eT>
976 inline
977 void
ones()978 diagview<eT>::ones()
979   {
980   arma_extra_debug_sigprint();
981 
982   (*this).fill(eT(1));
983   }
984 
985 
986 
987 template<typename eT>
988 inline
989 void
randu()990 diagview<eT>::randu()
991   {
992   arma_extra_debug_sigprint();
993 
994   Mat<eT>& x = const_cast< Mat<eT>& >(m);
995 
996   const uword local_n_elem = n_elem;
997 
998   for(uword ii=0; ii < local_n_elem; ++ii)
999     {
1000     x.at(ii+row_offset, ii+col_offset) = eT(arma_rng::randu<eT>());
1001     }
1002   }
1003 
1004 
1005 
1006 template<typename eT>
1007 inline
1008 void
randn()1009 diagview<eT>::randn()
1010   {
1011   arma_extra_debug_sigprint();
1012 
1013   Mat<eT>& x = const_cast< Mat<eT>& >(m);
1014 
1015   const uword local_n_elem = n_elem;
1016 
1017   for(uword ii=0; ii < local_n_elem; ++ii)
1018     {
1019     x.at(ii+row_offset, ii+col_offset) = eT(arma_rng::randn<eT>());
1020     }
1021   }
1022 
1023 
1024 
1025 //! @}
1026