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 SpSubview
20 //! @{
21 
22 
23 ///////////////////////////////////////////////////////////////////////////////
24 // SpSubview::iterator_base implementation                                   //
25 ///////////////////////////////////////////////////////////////////////////////
26 
27 template<typename eT>
28 inline
iterator_base(const SpSubview<eT> & in_M)29 SpSubview<eT>::iterator_base::iterator_base(const SpSubview<eT>& in_M)
30   : M(&in_M)
31   , internal_col(0)
32   , internal_pos(0)
33   {
34   // Technically this iterator is invalid (it may not point to a valid element)
35   }
36 
37 
38 
39 template<typename eT>
40 inline
iterator_base(const SpSubview<eT> & in_M,const uword in_col,const uword in_pos)41 SpSubview<eT>::iterator_base::iterator_base(const SpSubview<eT>& in_M, const uword in_col, const uword in_pos)
42   : M(&in_M)
43   , internal_col(in_col)
44   , internal_pos(in_pos)
45   {
46   // Nothing to do.
47   }
48 
49 
50 
51 ///////////////////////////////////////////////////////////////////////////////
52 // SpSubview::const_iterator implementation                                  //
53 ///////////////////////////////////////////////////////////////////////////////
54 
55 template<typename eT>
56 inline
const_iterator(const SpSubview<eT> & in_M,const uword initial_pos)57 SpSubview<eT>::const_iterator::const_iterator(const SpSubview<eT>& in_M, const uword initial_pos)
58   : iterator_base(in_M, 0, initial_pos)
59   {
60   // Corner case for empty subviews.
61   if(in_M.n_nonzero == 0)
62     {
63     iterator_base::internal_col = in_M.n_cols;
64     skip_pos                    = in_M.m.n_nonzero;
65     return;
66     }
67 
68   // Figure out the row and column of the position.
69   // lskip_pos holds the number of values which aren't part of this subview.
70   const uword aux_col = iterator_base::M->aux_col1;
71   const uword aux_row = iterator_base::M->aux_row1;
72   const uword ln_rows = iterator_base::M->n_rows;
73   const uword ln_cols = iterator_base::M->n_cols;
74 
75   uword cur_pos   = 0; // off by one because we might be searching for pos 0
76   uword lskip_pos = iterator_base::M->m.col_ptrs[aux_col];
77   uword cur_col   = 0;
78 
79   while(cur_pos < (iterator_base::internal_pos + 1))
80     {
81     // Have we stepped forward a column (or multiple columns)?
82     while(((lskip_pos + cur_pos) >= iterator_base::M->m.col_ptrs[cur_col + aux_col + 1]) && (cur_col < ln_cols))
83       {
84       ++cur_col;
85       }
86 
87     // See if the current position is in the subview.
88     const uword row_index = iterator_base::M->m.row_indices[cur_pos + lskip_pos];
89     if(row_index < aux_row)
90       {
91       ++lskip_pos; // not valid
92       }
93     else if(row_index < (aux_row + ln_rows))
94       {
95       ++cur_pos; // valid, in the subview
96       }
97     else
98       {
99       // skip to end of column
100       const uword next_colptr = iterator_base::M->m.col_ptrs[cur_col + aux_col + 1];
101       lskip_pos += (next_colptr - (cur_pos + lskip_pos));
102       }
103     }
104 
105   iterator_base::internal_col = cur_col;
106   skip_pos                    = lskip_pos;
107   }
108 
109 
110 
111 template<typename eT>
112 inline
const_iterator(const SpSubview<eT> & in_M,const uword in_row,const uword in_col)113 SpSubview<eT>::const_iterator::const_iterator(const SpSubview<eT>& in_M, const uword in_row, const uword in_col)
114   : iterator_base(in_M, in_col, 0)
115   {
116   // Corner case for empty subviews.
117   if(in_M.n_nonzero == 0)
118     {
119     // We must be at the last position.
120     iterator_base::internal_col = in_M.n_cols;
121     skip_pos                    = in_M.m.n_nonzero;
122     return;
123     }
124 
125   // We have a destination we want to be just after, but don't know what position that is.
126   // Because we have to count the points in this subview and not in this subview, this becomes a little difficult and slow.
127   const uword aux_col = iterator_base::M->aux_col1;
128   const uword aux_row = iterator_base::M->aux_row1;
129   const uword ln_rows = iterator_base::M->n_rows;
130   const uword ln_cols = iterator_base::M->n_cols;
131 
132   uword cur_pos = 0;
133   skip_pos = iterator_base::M->m.col_ptrs[aux_col];
134   uword cur_col = 0;
135 
136   // Skip any empty columns.
137   while(((skip_pos + cur_pos) >= iterator_base::M->m.col_ptrs[cur_col + aux_col + 1]) && (cur_col < ln_cols))
138     {
139     ++cur_col;
140     }
141 
142   while(cur_col < in_col)
143     {
144     // See if the current position is in the subview.
145     const uword row_index = iterator_base::M->m.row_indices[cur_pos + skip_pos];
146     if(row_index < aux_row)
147       {
148       ++skip_pos;
149       }
150     else if(row_index < (aux_row + ln_rows))
151       {
152       ++cur_pos;
153       }
154     else
155       {
156       // skip to end of column
157       const uword next_colptr = iterator_base::M->m.col_ptrs[cur_col + aux_col + 1];
158       skip_pos += (next_colptr - (cur_pos + skip_pos));
159       }
160 
161     // Have we stepped forward a column (or multiple columns)?
162     while(((skip_pos + cur_pos) >= iterator_base::M->m.col_ptrs[cur_col + aux_col + 1]) && (cur_col < ln_cols))
163       {
164       ++cur_col;
165       }
166     }
167 
168   // Now we are either on the right column or ahead of it.
169   if(cur_col == in_col)
170     {
171     // We have to find the right row index.
172     uword row_index = iterator_base::M->m.row_indices[cur_pos + skip_pos];
173     while((row_index < (in_row + aux_row)))
174       {
175       if(row_index < aux_row)
176         {
177         ++skip_pos;
178         }
179       else
180         {
181         ++cur_pos;
182         }
183 
184       // Ensure we didn't step forward a column; if we did, we need to stop.
185       while(((skip_pos + cur_pos) >= iterator_base::M->m.col_ptrs[cur_col + aux_col + 1]) && (cur_col < ln_cols))
186         {
187         ++cur_col;
188         }
189 
190       if(cur_col != in_col)
191         {
192         break;
193         }
194 
195       row_index = iterator_base::M->m.row_indices[cur_pos + skip_pos];
196       }
197     }
198 
199   // Now we need to find the next valid position in the subview.
200   uword row_index;
201   while(true)
202     {
203     const uword next_colptr = iterator_base::M->m.col_ptrs[cur_col + aux_col + 1];
204     row_index = iterator_base::M->m.row_indices[cur_pos + skip_pos];
205 
206     // Are we at the last position?
207     if(cur_col >= ln_cols)
208       {
209       cur_col = ln_cols;
210       // Make sure we will be pointing at the last element in the parent matrix.
211       skip_pos = iterator_base::M->m.n_nonzero - iterator_base::M->n_nonzero;
212       break;
213       }
214 
215     if(row_index < aux_row)
216       {
217       ++skip_pos;
218       }
219     else if(row_index < (aux_row + ln_rows))
220       {
221       break; // found
222       }
223     else
224       {
225       skip_pos += (next_colptr - (cur_pos + skip_pos));
226       }
227 
228     // Did we move any columns?
229     while(((skip_pos + cur_pos) >= iterator_base::M->m.col_ptrs[cur_col + aux_col + 1]) && (cur_col < ln_cols))
230       {
231       ++cur_col;
232       }
233     }
234 
235   // It is possible we have moved another column.
236   while(((skip_pos + cur_pos) >= iterator_base::M->m.col_ptrs[cur_col + aux_col + 1]) && (cur_col < ln_cols))
237     {
238     ++cur_col;
239     }
240 
241   iterator_base::internal_pos = cur_pos;
242   iterator_base::internal_col = cur_col;
243   }
244 
245 
246 
247 template<typename eT>
248 inline
const_iterator(const SpSubview<eT> & in_M,uword in_row,uword in_col,uword in_pos,uword in_skip_pos)249 SpSubview<eT>::const_iterator::const_iterator(const SpSubview<eT>& in_M, uword in_row, uword in_col, uword in_pos, uword in_skip_pos)
250   : iterator_base(in_M, in_col, in_pos)
251   , skip_pos(in_skip_pos)
252   {
253   arma_ignore(in_row);
254 
255   // Nothing to do.
256   }
257 
258 
259 
260 template<typename eT>
261 inline
const_iterator(const const_iterator & other)262 SpSubview<eT>::const_iterator::const_iterator(const const_iterator& other)
263   : iterator_base(*other.M, other.internal_col, other.internal_pos)
264   , skip_pos(other.skip_pos)
265   {
266   // Nothing to do.
267   }
268 
269 
270 
271 template<typename eT>
272 arma_inline
273 eT
operator *() const274 SpSubview<eT>::const_iterator::operator*() const
275   {
276   return iterator_base::M->m.values[iterator_base::internal_pos + skip_pos];
277   }
278 
279 
280 
281 template<typename eT>
282 inline
283 arma_hot
284 typename SpSubview<eT>::const_iterator&
operator ++()285 SpSubview<eT>::const_iterator::operator++()
286   {
287   const uword aux_col = iterator_base::M->aux_col1;
288   const uword aux_row = iterator_base::M->aux_row1;
289   const uword ln_rows = iterator_base::M->n_rows;
290   const uword ln_cols = iterator_base::M->n_cols;
291 
292   uword cur_col   = iterator_base::internal_col;
293   uword cur_pos   = iterator_base::internal_pos + 1;
294   uword lskip_pos = skip_pos;
295   uword row_index;
296 
297   while(true)
298     {
299     const uword next_colptr = iterator_base::M->m.col_ptrs[cur_col + aux_col + 1];
300     row_index = iterator_base::M->m.row_indices[cur_pos + lskip_pos];
301 
302     // Did we move any columns?
303     while((cur_col < ln_cols) && ((lskip_pos + cur_pos) >= iterator_base::M->m.col_ptrs[cur_col + aux_col + 1]))
304       {
305       ++cur_col;
306       }
307 
308     // Are we at the last position?
309     if(cur_col >= ln_cols)
310       {
311       cur_col = ln_cols;
312       // Make sure we will be pointing at the last element in the parent matrix.
313       lskip_pos = iterator_base::M->m.n_nonzero - iterator_base::M->n_nonzero;
314       break;
315       }
316 
317     if(row_index < aux_row)
318       {
319       ++lskip_pos;
320       }
321     else if(row_index < (aux_row + ln_rows))
322       {
323       break; // found
324       }
325     else
326       {
327       lskip_pos += (next_colptr - (cur_pos + lskip_pos));
328       }
329     }
330 
331   iterator_base::internal_pos = cur_pos;
332   iterator_base::internal_col = cur_col;
333   skip_pos                    = lskip_pos;
334 
335   return *this;
336   }
337 
338 
339 
340 template<typename eT>
341 inline
342 arma_warn_unused
343 typename SpSubview<eT>::const_iterator
operator ++(int)344 SpSubview<eT>::const_iterator::operator++(int)
345   {
346   typename SpSubview<eT>::const_iterator tmp(*this);
347 
348   ++(*this);
349 
350   return tmp;
351   }
352 
353 
354 
355 template<typename eT>
356 inline
357 arma_hot
358 typename SpSubview<eT>::const_iterator&
operator --()359 SpSubview<eT>::const_iterator::operator--()
360   {
361   const uword aux_col = iterator_base::M->aux_col1;
362   const uword aux_row = iterator_base::M->aux_row1;
363   const uword ln_rows = iterator_base::M->n_rows;
364 
365   uword cur_col  = iterator_base::internal_col;
366   uword cur_pos  = iterator_base::internal_pos - 1;
367 
368   // Special condition for end of iterator.
369   if((skip_pos + cur_pos + 1) == iterator_base::M->m.n_nonzero)
370     {
371     // We are at the last element.  So we need to set skip_pos back to what it
372     // would be if we didn't manually modify it back in operator++().
373     skip_pos = iterator_base::M->m.col_ptrs[cur_col + aux_col] - iterator_base::internal_pos;
374     }
375 
376   uword row_index;
377 
378   while(true)
379     {
380     const uword colptr = iterator_base::M->m.col_ptrs[cur_col + aux_col];
381     row_index = iterator_base::M->m.row_indices[cur_pos + skip_pos];
382 
383     // Did we move back any columns?
384     while((skip_pos + cur_pos) < iterator_base::M->m.col_ptrs[cur_col + aux_col])
385       {
386       --cur_col;
387       }
388 
389     if(row_index < aux_row)
390       {
391       skip_pos -= (colptr - (cur_pos + skip_pos) + 1);
392       }
393     else if(row_index < (aux_row + ln_rows))
394       {
395       break; // found
396       }
397     else
398       {
399       --skip_pos;
400       }
401     }
402 
403   iterator_base::internal_pos = cur_pos;
404   iterator_base::internal_col = cur_col;
405 
406   return *this;
407   }
408 
409 
410 
411 template<typename eT>
412 inline
413 arma_warn_unused
414 typename SpSubview<eT>::const_iterator
operator --(int)415 SpSubview<eT>::const_iterator::operator--(int)
416   {
417   typename SpSubview<eT>::const_iterator tmp(*this);
418 
419   --(*this);
420 
421   return tmp;
422   }
423 
424 
425 
426 template<typename eT>
427 inline
428 arma_hot
429 bool
operator ==(const const_iterator & rhs) const430 SpSubview<eT>::const_iterator::operator==(const const_iterator& rhs) const
431   {
432   return (rhs.row() == (*this).row()) && (rhs.col() == iterator_base::internal_col);
433   }
434 
435 
436 
437 template<typename eT>
438 inline
439 arma_hot
440 bool
operator !=(const const_iterator & rhs) const441 SpSubview<eT>::const_iterator::operator!=(const const_iterator& rhs) const
442   {
443   return (rhs.row() != (*this).row()) || (rhs.col() != iterator_base::internal_col);
444   }
445 
446 
447 
448 template<typename eT>
449 inline
450 arma_hot
451 bool
operator ==(const typename SpMat<eT>::const_iterator & rhs) const452 SpSubview<eT>::const_iterator::operator==(const typename SpMat<eT>::const_iterator& rhs) const
453   {
454   return (rhs.row() == (*this).row()) && (rhs.col() == iterator_base::internal_col);
455   }
456 
457 
458 
459 template<typename eT>
460 inline
461 arma_hot
462 bool
operator !=(const typename SpMat<eT>::const_iterator & rhs) const463 SpSubview<eT>::const_iterator::operator!=(const typename SpMat<eT>::const_iterator& rhs) const
464   {
465   return (rhs.row() != (*this).row()) || (rhs.col() != iterator_base::internal_col);
466   }
467 
468 
469 
470 template<typename eT>
471 inline
472 arma_hot
473 bool
operator ==(const const_row_iterator & rhs) const474 SpSubview<eT>::const_iterator::operator==(const const_row_iterator& rhs) const
475   {
476   return (rhs.row() == (*this).row()) && (rhs.col() == iterator_base::internal_col);
477   }
478 
479 
480 
481 template<typename eT>
482 inline
483 arma_hot
484 bool
operator !=(const const_row_iterator & rhs) const485 SpSubview<eT>::const_iterator::operator!=(const const_row_iterator& rhs) const
486   {
487   return (rhs.row() != (*this).row()) || (rhs.col() != iterator_base::internal_col);
488   }
489 
490 
491 
492 template<typename eT>
493 inline
494 arma_hot
495 bool
operator ==(const typename SpMat<eT>::const_row_iterator & rhs) const496 SpSubview<eT>::const_iterator::operator==(const typename SpMat<eT>::const_row_iterator& rhs) const
497   {
498   return (rhs.row() == (*this).row()) && (rhs.col() == iterator_base::internal_col);
499   }
500 
501 
502 
503 template<typename eT>
504 inline
505 arma_hot
506 bool
operator !=(const typename SpMat<eT>::const_row_iterator & rhs) const507 SpSubview<eT>::const_iterator::operator!=(const typename SpMat<eT>::const_row_iterator& rhs) const
508   {
509   return (rhs.row() != (*this).row()) || (rhs.col() != iterator_base::internal_col);
510   }
511 
512 
513 
514 ///////////////////////////////////////////////////////////////////////////////
515 // SpSubview<eT>::iterator implementation                                    //
516 ///////////////////////////////////////////////////////////////////////////////
517 
518 template<typename eT>
519 inline
520 arma_hot
521 SpValProxy< SpSubview<eT> >
operator *()522 SpSubview<eT>::iterator::operator*()
523   {
524   return SpValProxy< SpSubview<eT> >(
525     const_iterator::row(),
526     iterator_base::col(),
527     access::rw(*iterator_base::M),
528     &(access::rw(iterator_base::M->m.values[iterator_base::internal_pos + const_iterator::skip_pos])));
529   }
530 
531 
532 
533 template<typename eT>
534 inline
535 arma_hot
536 typename SpSubview<eT>::iterator&
operator ++()537 SpSubview<eT>::iterator::operator++()
538   {
539   const_iterator::operator++();
540   return *this;
541   }
542 
543 
544 
545 template<typename eT>
546 inline
547 arma_warn_unused
548 typename SpSubview<eT>::iterator
operator ++(int)549 SpSubview<eT>::iterator::operator++(int)
550   {
551   typename SpSubview<eT>::iterator tmp(*this);
552 
553   const_iterator::operator++();
554 
555   return tmp;
556   }
557 
558 
559 
560 template<typename eT>
561 inline
562 arma_hot
563 typename SpSubview<eT>::iterator&
operator --()564 SpSubview<eT>::iterator::operator--()
565   {
566   const_iterator::operator--();
567   return *this;
568   }
569 
570 
571 
572 template<typename eT>
573 inline
574 arma_warn_unused
575 typename SpSubview<eT>::iterator
operator --(int)576 SpSubview<eT>::iterator::operator--(int)
577   {
578   typename SpSubview<eT>::iterator tmp(*this);
579 
580   const_iterator::operator--();
581 
582   return tmp;
583   }
584 
585 
586 
587 ///////////////////////////////////////////////////////////////////////////////
588 // SpSubview<eT>::const_row_iterator implementation                          //
589 ///////////////////////////////////////////////////////////////////////////////
590 
591 template<typename eT>
592 inline
const_row_iterator()593 SpSubview<eT>::const_row_iterator::const_row_iterator()
594   : iterator_base()
595   , internal_row(0)
596   , actual_pos(0)
597   {
598   }
599 
600 
601 
602 template<typename eT>
603 inline
const_row_iterator(const SpSubview<eT> & in_M,uword initial_pos)604 SpSubview<eT>::const_row_iterator::const_row_iterator(const SpSubview<eT>& in_M, uword initial_pos)
605   : iterator_base(in_M, 0, initial_pos)
606   , internal_row(0)
607   , actual_pos(0)
608   {
609   // Corner case for the end of a subview.
610   if(initial_pos == in_M.n_nonzero)
611     {
612     iterator_base::internal_col = 0;
613     internal_row = in_M.n_rows;
614     return;
615     }
616 
617   const uword aux_col = iterator_base::M->aux_col1;
618   const uword aux_row = iterator_base::M->aux_row1;
619 
620   // We don't count zeros in our position count, so we have to find the nonzero
621   // value corresponding to the given initial position, and we also have to skip
622   // any nonzero elements that aren't a part of the subview.
623 
624   uword cur_pos = std::numeric_limits<uword>::max();
625   uword cur_actual_pos = 0;
626 
627   // Since we don't know where the elements are in each row, we have to loop
628   // across all columns looking for elements in row 0 and add to our sum, then
629   // in row 1, and so forth, until we get to the desired position.
630   for(uword row = 0; row < iterator_base::M->n_rows; ++row)
631     {
632     for(uword col = 0; col < iterator_base::M->n_cols; ++col)
633       {
634       // Find the first element with row greater than or equal to row + aux_row.
635       const uword      col_offset = iterator_base::M->m.col_ptrs[col + aux_col    ];
636       const uword next_col_offset = iterator_base::M->m.col_ptrs[col + aux_col + 1];
637 
638       const uword* start_ptr = &iterator_base::M->m.row_indices[     col_offset];
639       const uword*   end_ptr = &iterator_base::M->m.row_indices[next_col_offset];
640 
641       if(start_ptr != end_ptr)
642         {
643         const uword* pos_ptr = std::lower_bound(start_ptr, end_ptr, row + aux_row);
644 
645         const uword offset = uword(pos_ptr - start_ptr);
646 
647         if(iterator_base::M->m.row_indices[col_offset + offset] == row + aux_row)
648           {
649           cur_actual_pos = col_offset + offset;
650 
651           // Increment position portably.
652           if(cur_pos == std::numeric_limits<uword>::max())
653             cur_pos = 0;
654           else
655             ++cur_pos;
656 
657           // Do we terminate?
658           if(cur_pos == initial_pos)
659             {
660             internal_row = row;
661             iterator_base::internal_col = col;
662             iterator_base::internal_pos = cur_pos;
663             actual_pos = cur_actual_pos;
664 
665             return;
666             }
667           }
668         }
669       }
670     }
671 
672   // This shouldn't happen.
673   iterator_base::internal_pos = iterator_base::M->n_nonzero;
674   iterator_base::internal_col = 0;
675   internal_row = iterator_base::M->n_rows;
676   actual_pos = iterator_base::M->n_nonzero;
677   }
678 
679 
680 
681 template<typename eT>
682 inline
const_row_iterator(const SpSubview<eT> & in_M,uword in_row,uword in_col)683 SpSubview<eT>::const_row_iterator::const_row_iterator(const SpSubview<eT>& in_M, uword in_row, uword in_col)
684   : iterator_base(in_M, in_col, 0)
685   , internal_row(0)
686   , actual_pos(0)
687   {
688   // Start our search in the given row.  We need to find two things:
689   //
690   //  1. The first nonzero element (iterating by rows) after (in_row, in_col).
691   //  2. The number of nonzero elements (iterating by rows) that come before
692   //     (in_row, in_col).
693   //
694   // We'll find these simultaneously, though we will have to loop over all
695   // columns.
696 
697   const uword aux_col = iterator_base::M->aux_col1;
698   const uword aux_row = iterator_base::M->aux_row1;
699 
700   // This will hold the total number of points in the subview with rows less
701   // than in_row.
702   uword cur_pos = 0;
703   uword cur_min_row = iterator_base::M->n_rows;
704   uword cur_min_col = 0;
705   uword cur_actual_pos = 0;
706 
707   for(uword col = 0; col < iterator_base::M->n_cols; ++col)
708     {
709     // Find the first element with row greater than or equal to in_row.
710     const uword      col_offset = iterator_base::M->m.col_ptrs[col + aux_col    ];
711     const uword next_col_offset = iterator_base::M->m.col_ptrs[col + aux_col + 1];
712 
713     const uword* start_ptr = &iterator_base::M->m.row_indices[     col_offset];
714     const uword*   end_ptr = &iterator_base::M->m.row_indices[next_col_offset];
715 
716     if(start_ptr != end_ptr)
717       {
718       // First let us find the first element that is in the subview.
719       const uword* first_subview_ptr = std::lower_bound(start_ptr, end_ptr, aux_row);
720 
721       if(first_subview_ptr != end_ptr && (*first_subview_ptr) < aux_row + iterator_base::M->n_rows)
722         {
723         // There exists at least one element in the subview.
724         const uword* pos_ptr = std::lower_bound(first_subview_ptr, end_ptr, aux_row + in_row);
725 
726         // This is the number of elements in the subview with row index less
727         // than in_row.
728         cur_pos += uword(pos_ptr - first_subview_ptr);
729 
730         if(pos_ptr != end_ptr && (*pos_ptr) < aux_row + iterator_base::M->n_rows)
731           {
732           // This is the row index of the first element in the column with row
733           // index greater than or equal to in_row + aux_row.
734           if((*pos_ptr) - aux_row < cur_min_row)
735             {
736             // If we are in the desired row but before the desired column, we
737             // can't take this.
738             if(col >= in_col)
739               {
740               cur_min_row = (*pos_ptr) - aux_row;
741               cur_min_col = col;
742               cur_actual_pos = col_offset + (pos_ptr - start_ptr);
743               }
744             }
745           }
746         }
747       }
748     }
749 
750   // Now we know what the minimum row is.
751   internal_row = cur_min_row;
752   iterator_base::internal_col = cur_min_col;
753   iterator_base::internal_pos = cur_pos;
754   actual_pos = cur_actual_pos;
755   }
756 
757 
758 
759 template<typename eT>
760 inline
const_row_iterator(const const_row_iterator & other)761 SpSubview<eT>::const_row_iterator::const_row_iterator(const const_row_iterator& other)
762   : iterator_base(*other.M, other.internal_col, other.internal_pos)
763   , internal_row(other.internal_row)
764   , actual_pos(other.actual_pos)
765   {
766   // Nothing to do.
767   }
768 
769 
770 
771 template<typename eT>
772 inline
773 arma_hot
774 typename SpSubview<eT>::const_row_iterator&
operator ++()775 SpSubview<eT>::const_row_iterator::operator++()
776   {
777   // We just need to find the next nonzero element.
778   ++iterator_base::internal_pos;
779 
780   // If we have exceeded the bounds, update accordingly.
781   if(iterator_base::internal_pos >= iterator_base::M->n_nonzero)
782     {
783     internal_row = iterator_base::M->n_rows;
784     iterator_base::internal_col = 0;
785     actual_pos = iterator_base::M->n_nonzero;
786 
787     return *this;
788     }
789 
790   const uword aux_col  = iterator_base::M->aux_col1;
791   const uword aux_row  = iterator_base::M->aux_row1;
792   const uword M_n_cols = iterator_base::M->n_cols;
793 
794   // Otherwise, we need to search.  We have to loop over all of the columns in
795   // the subview.
796   uword next_min_row = iterator_base::M->n_rows;
797   uword next_min_col = 0;
798   uword next_actual_pos = 0;
799 
800   for(uword col = iterator_base::internal_col + 1; col < M_n_cols; ++col)
801     {
802     // Find the first element with row greater than or equal to row.
803     const uword      col_offset = iterator_base::M->m.col_ptrs[col + aux_col    ];
804     const uword next_col_offset = iterator_base::M->m.col_ptrs[col + aux_col + 1];
805 
806     const uword* start_ptr = &iterator_base::M->m.row_indices[     col_offset];
807     const uword* end_ptr   = &iterator_base::M->m.row_indices[next_col_offset];
808 
809     if(start_ptr != end_ptr)
810       {
811       // Find the first element in the column with row greater than or equal to
812       // the current row.  Since this is a subview, it's possible that we may
813       // find rows past the end of the subview.
814       const uword* pos_ptr = std::lower_bound(start_ptr, end_ptr, internal_row + aux_row);
815 
816       if(pos_ptr != end_ptr)
817         {
818         // We found something; is the row index correct?
819         if((*pos_ptr) == internal_row + aux_row && (*pos_ptr) < aux_row + iterator_base::M->n_rows)
820           {
821           // Exact match---so we are done.
822           iterator_base::internal_col = col;
823           actual_pos = col_offset + (pos_ptr - start_ptr);
824           return *this;
825           }
826         else if((*pos_ptr) < next_min_row + aux_row && (*pos_ptr) < aux_row + iterator_base::M->n_rows)
827           {
828           // The first element in this column is in a subsequent row, but it's
829           // the minimum row we've seen so far.
830           next_min_row = (*pos_ptr) - aux_row;
831           next_min_col = col;
832           next_actual_pos = col_offset + (pos_ptr - start_ptr);
833           }
834         else if((*pos_ptr) == next_min_row + aux_row && col < next_min_col && (*pos_ptr) < aux_row + iterator_base::M->n_rows)
835           {
836           // The first element in this column is in a subsequent row that we
837           // already have another elemnt for, but the column index is less so
838           // this element will come first.
839           next_min_col = col;
840           next_actual_pos = col_offset + (pos_ptr - start_ptr);
841           }
842         }
843       }
844     }
845 
846   // Restart the search in the next row.
847   for(uword col = 0; col <= iterator_base::internal_col; ++col)
848     {
849     // Find the first element with row greater than or equal to row + 1.
850     const uword      col_offset = iterator_base::M->m.col_ptrs[col + aux_col    ];
851     const uword next_col_offset = iterator_base::M->m.col_ptrs[col + aux_col + 1];
852 
853     const uword* start_ptr = &iterator_base::M->m.row_indices[     col_offset];
854     const uword* end_ptr   = &iterator_base::M->m.row_indices[next_col_offset];
855 
856     if(start_ptr != end_ptr)
857       {
858       const uword* pos_ptr = std::lower_bound(start_ptr, end_ptr, internal_row + aux_row + 1);
859 
860       if(pos_ptr != end_ptr)
861         {
862         // We found something in the column, but is the row index correct?
863         if((*pos_ptr) == internal_row + aux_row + 1 && (*pos_ptr) < aux_row + iterator_base::M->n_rows)
864           {
865           // Exact match---so we are done.
866           iterator_base::internal_col = col;
867           internal_row++;
868           actual_pos = col_offset + (pos_ptr - start_ptr);
869           return *this;
870           }
871         else if((*pos_ptr) < next_min_row + aux_row && (*pos_ptr) < aux_row + iterator_base::M->n_rows)
872           {
873           // The first element in this column is in a subsequent row, but it's
874           // the minimum row we've seen so far.
875           next_min_row = (*pos_ptr) - aux_row;
876           next_min_col = col;
877           next_actual_pos = col_offset + (pos_ptr - start_ptr);
878           }
879         else if((*pos_ptr) == next_min_row + aux_row && col < next_min_col && (*pos_ptr) < aux_row + iterator_base::M->n_rows)
880           {
881           // We've found a better column.
882           next_min_col = col;
883           next_actual_pos = col_offset + (pos_ptr - start_ptr);
884           }
885         }
886       }
887     }
888 
889   iterator_base::internal_col = next_min_col;
890   internal_row = next_min_row;
891   actual_pos = next_actual_pos;
892 
893   return *this;
894   }
895 
896 
897 
898 template<typename eT>
899 inline
900 arma_warn_unused
901 typename SpSubview<eT>::const_row_iterator
operator ++(int)902 SpSubview<eT>::const_row_iterator::operator++(int)
903   {
904   typename SpSubview<eT>::const_row_iterator tmp(*this);
905 
906   ++(*this);
907 
908   return tmp;
909   }
910 
911 
912 
913 template<typename eT>
914 inline
915 arma_hot
916 typename SpSubview<eT>::const_row_iterator&
operator --()917 SpSubview<eT>::const_row_iterator::operator--()
918   {
919   if(iterator_base::internal_pos == 0)
920     {
921     // We are already at the beginning.
922     return *this;
923     }
924 
925   iterator_base::internal_pos--;
926 
927   const uword aux_col = iterator_base::M->aux_col1;
928   const uword aux_row = iterator_base::M->aux_row1;
929 
930   // We have to search backwards.
931   uword max_row = 0;
932   uword max_col = 0;
933   uword next_actual_pos = 0;
934 
935   for(uword col = iterator_base::internal_col; col >= 1; --col)
936     {
937     // Find the first element with row greater than or equal to in_row + 1.
938     const uword      col_offset = iterator_base::M->m.col_ptrs[col + aux_col - 1];
939     const uword next_col_offset = iterator_base::M->m.col_ptrs[col + aux_col   ];
940 
941     const uword* start_ptr = &iterator_base::M->m.row_indices[     col_offset];
942     const uword* end_ptr   = &iterator_base::M->m.row_indices[next_col_offset];
943 
944     if(start_ptr != end_ptr)
945       {
946       // There are elements in this column.
947       const uword* pos_ptr = std::lower_bound(start_ptr, end_ptr, internal_row + aux_row + 1);
948 
949       if(pos_ptr != start_ptr)
950         {
951         if(*(pos_ptr - 1) > max_row + aux_row)
952           {
953           // There are elements in this column with row index < internal_row.
954           max_row = *(pos_ptr - 1) - aux_row;
955           max_col = col - 1;
956           next_actual_pos = col_offset + (pos_ptr - 1 - start_ptr);
957           }
958         else if(*(pos_ptr - 1) == max_row + aux_row && (col - 1) >= max_col)
959           {
960           max_col = col - 1;
961           next_actual_pos = col_offset + (pos_ptr - 1 - start_ptr);
962           }
963         }
964       }
965     }
966 
967   for(uword col = iterator_base::M->n_cols - 1; col >= iterator_base::internal_col; --col)
968     {
969     // Find the first element with row greater than or equal to row + 1.
970     const uword      col_offset = iterator_base::M->m.col_ptrs[col + aux_col    ];
971     const uword next_col_offset = iterator_base::M->m.col_ptrs[col + aux_col + 1];
972 
973     const uword* start_ptr = &iterator_base::M->m.row_indices[     col_offset];
974     const uword*   end_ptr = &iterator_base::M->m.row_indices[next_col_offset];
975 
976     if(start_ptr != end_ptr)
977       {
978       // There are elements in this column.
979       const uword* pos_ptr = std::lower_bound(start_ptr, end_ptr, internal_row + aux_row);
980 
981       if(pos_ptr != start_ptr)
982         {
983         // There are elements in this column with row index < internal_row.
984         if(*(pos_ptr - 1) > max_row + aux_row)
985           {
986           max_row = *(pos_ptr - 1) - aux_row;
987           max_col = col;
988           next_actual_pos = col_offset + (pos_ptr - 1 - start_ptr);
989           }
990         else if(*(pos_ptr - 1) == max_row + aux_row && col >= max_col)
991           {
992           max_col = col;
993           next_actual_pos = col_offset + (pos_ptr - 1 - start_ptr);
994           }
995         }
996       }
997 
998     if(col == 0) // Catch edge case that the loop termination condition won't.
999       {
1000       break;
1001       }
1002     }
1003 
1004   iterator_base::internal_col = max_col;
1005   internal_row = max_row;
1006   actual_pos = next_actual_pos;
1007 
1008   return *this;
1009   }
1010 
1011 
1012 
1013 template<typename eT>
1014 inline
1015 arma_warn_unused
1016 typename SpSubview<eT>::const_row_iterator
operator --(int)1017 SpSubview<eT>::const_row_iterator::operator--(int)
1018   {
1019   typename SpSubview<eT>::const_row_iterator tmp(*this);
1020 
1021   --(*this);
1022 
1023   return tmp;
1024   }
1025 
1026 
1027 
1028 template<typename eT>
1029 inline
1030 arma_hot
1031 bool
operator ==(const const_iterator & rhs) const1032 SpSubview<eT>::const_row_iterator::operator==(const const_iterator& rhs) const
1033   {
1034   return (rhs.row() == row()) && (rhs.col() == iterator_base::internal_col);
1035   }
1036 
1037 
1038 
1039 template<typename eT>
1040 inline
1041 arma_hot
1042 bool
operator !=(const const_iterator & rhs) const1043 SpSubview<eT>::const_row_iterator::operator!=(const const_iterator& rhs) const
1044   {
1045   return (rhs.row() != row()) || (rhs.col() != iterator_base::internal_col);
1046   }
1047 
1048 
1049 
1050 template<typename eT>
1051 inline
1052 arma_hot
1053 bool
operator ==(const typename SpMat<eT>::const_iterator & rhs) const1054 SpSubview<eT>::const_row_iterator::operator==(const typename SpMat<eT>::const_iterator& rhs) const
1055   {
1056   return (rhs.row() == row()) && (rhs.col() == iterator_base::internal_col);
1057   }
1058 
1059 
1060 
1061 template<typename eT>
1062 inline
1063 arma_hot
1064 bool
operator !=(const typename SpMat<eT>::const_iterator & rhs) const1065 SpSubview<eT>::const_row_iterator::operator!=(const typename SpMat<eT>::const_iterator& rhs) const
1066   {
1067   return (rhs.row() != row()) || (rhs.col() != iterator_base::internal_col);
1068   }
1069 
1070 
1071 
1072 template<typename eT>
1073 inline
1074 arma_hot
1075 bool
operator ==(const const_row_iterator & rhs) const1076 SpSubview<eT>::const_row_iterator::operator==(const const_row_iterator& rhs) const
1077   {
1078   return (rhs.row() == row()) && (rhs.col() == iterator_base::internal_col);
1079   }
1080 
1081 
1082 
1083 template<typename eT>
1084 inline
1085 arma_hot
1086 bool
operator !=(const const_row_iterator & rhs) const1087 SpSubview<eT>::const_row_iterator::operator!=(const const_row_iterator& rhs) const
1088   {
1089   return (rhs.row() != row()) || (rhs.col() != iterator_base::internal_col);
1090   }
1091 
1092 
1093 
1094 template<typename eT>
1095 inline
1096 arma_hot
1097 bool
operator ==(const typename SpMat<eT>::const_row_iterator & rhs) const1098 SpSubview<eT>::const_row_iterator::operator==(const typename SpMat<eT>::const_row_iterator& rhs) const
1099   {
1100   return (rhs.row() == row()) && (rhs.col() == iterator_base::internal_col);
1101   }
1102 
1103 
1104 
1105 template<typename eT>
1106 inline
1107 arma_hot
1108 bool
operator !=(const typename SpMat<eT>::const_row_iterator & rhs) const1109 SpSubview<eT>::const_row_iterator::operator!=(const typename SpMat<eT>::const_row_iterator& rhs) const
1110   {
1111   return (rhs.row() != row()) || (rhs.col() != iterator_base::internal_col);
1112   }
1113 
1114 
1115 
1116 ///////////////////////////////////////////////////////////////////////////////
1117 // SpSubview<eT>::row_iterator implementation                                //
1118 ///////////////////////////////////////////////////////////////////////////////
1119 
1120 template<typename eT>
1121 inline
1122 arma_hot
1123 SpValProxy< SpSubview<eT> >
operator *()1124 SpSubview<eT>::row_iterator::operator*()
1125   {
1126   return SpValProxy< SpSubview<eT> >(
1127     const_row_iterator::internal_row,
1128     iterator_base::internal_col,
1129     access::rw(*iterator_base::M),
1130     &access::rw(iterator_base::M->m.values[const_row_iterator::actual_pos]));
1131   }
1132 
1133 
1134 
1135 template<typename eT>
1136 inline
1137 arma_hot
1138 typename SpSubview<eT>::row_iterator&
operator ++()1139 SpSubview<eT>::row_iterator::operator++()
1140   {
1141   const_row_iterator::operator++();
1142   return *this;
1143   }
1144 
1145 
1146 
1147 template<typename eT>
1148 inline
1149 arma_warn_unused
1150 typename SpSubview<eT>::row_iterator
operator ++(int)1151 SpSubview<eT>::row_iterator::operator++(int)
1152   {
1153   typename SpSubview<eT>::row_iterator tmp(*this);
1154 
1155   ++(*this);
1156 
1157   return tmp;
1158   }
1159 
1160 
1161 
1162 template<typename eT>
1163 inline
1164 arma_hot
1165 typename SpSubview<eT>::row_iterator&
operator --()1166 SpSubview<eT>::row_iterator::operator--()
1167   {
1168   const_row_iterator::operator--();
1169   return *this;
1170   }
1171 
1172 
1173 
1174 template<typename eT>
1175 inline
1176 arma_warn_unused
1177 typename SpSubview<eT>::row_iterator
operator --(int)1178 SpSubview<eT>::row_iterator::operator--(int)
1179   {
1180   typename SpSubview<eT>::row_iterator tmp(*this);
1181 
1182   --(*this);
1183 
1184   return tmp;
1185   }
1186 
1187 
1188 //! @}
1189