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 arrayops
20 //! @{
21
22
23
24 template<typename eT>
25 arma_inline
26 void
copy(eT * dest,const eT * src,const uword n_elem)27 arrayops::copy(eT* dest, const eT* src, const uword n_elem)
28 {
29 if( (dest == src) || (n_elem == 0) ) { return; }
30
31 if(is_cx<eT>::no)
32 {
33 if(n_elem <= 9)
34 {
35 arrayops::copy_small(dest, src, n_elem);
36 }
37 else
38 {
39 std::memcpy(dest, src, n_elem*sizeof(eT));
40 }
41 }
42 else
43 {
44 std::memcpy(dest, src, n_elem*sizeof(eT));
45 }
46 }
47
48
49
50 template<typename eT>
51 arma_cold
52 inline
53 void
copy_small(eT * dest,const eT * src,const uword n_elem)54 arrayops::copy_small(eT* dest, const eT* src, const uword n_elem)
55 {
56 switch(n_elem)
57 {
58 case 9: dest[ 8] = src[ 8];
59 // fallthrough
60 case 8: dest[ 7] = src[ 7];
61 // fallthrough
62 case 7: dest[ 6] = src[ 6];
63 // fallthrough
64 case 6: dest[ 5] = src[ 5];
65 // fallthrough
66 case 5: dest[ 4] = src[ 4];
67 // fallthrough
68 case 4: dest[ 3] = src[ 3];
69 // fallthrough
70 case 3: dest[ 2] = src[ 2];
71 // fallthrough
72 case 2: dest[ 1] = src[ 1];
73 // fallthrough
74 case 1: dest[ 0] = src[ 0];
75 // fallthrough
76 default: ;
77 }
78 }
79
80
81
82 template<typename eT>
83 inline
84 void
fill_zeros(eT * dest,const uword n_elem)85 arrayops::fill_zeros(eT* dest, const uword n_elem)
86 {
87 typedef typename get_pod_type<eT>::result pod_type;
88
89 if(n_elem == 0) { return; }
90
91 if(std::numeric_limits<eT>::is_integer || std::numeric_limits<pod_type>::is_iec559)
92 {
93 std::memset((void*)dest, 0, sizeof(eT)*n_elem);
94 }
95 else
96 {
97 arrayops::inplace_set_simple(dest, eT(0), n_elem);
98 }
99 }
100
101
102
103 template<typename eT>
104 arma_hot
105 inline
106 void
replace(eT * mem,const uword n_elem,const eT old_val,const eT new_val)107 arrayops::replace(eT* mem, const uword n_elem, const eT old_val, const eT new_val)
108 {
109 if(arma_isnan(old_val))
110 {
111 for(uword i=0; i<n_elem; ++i)
112 {
113 eT& val = mem[i];
114
115 val = (arma_isnan(val)) ? new_val : val;
116 }
117 }
118 else
119 {
120 for(uword i=0; i<n_elem; ++i)
121 {
122 eT& val = mem[i];
123
124 val = (val == old_val) ? new_val : val;
125 }
126 }
127 }
128
129
130
131 template<typename eT>
132 arma_hot
133 inline
134 void
clean(eT * mem,const uword n_elem,const eT abs_limit,const typename arma_not_cx<eT>::result * junk)135 arrayops::clean(eT* mem, const uword n_elem, const eT abs_limit, const typename arma_not_cx<eT>::result* junk)
136 {
137 arma_ignore(junk);
138
139 for(uword i=0; i<n_elem; ++i)
140 {
141 eT& val = mem[i];
142
143 val = (eop_aux::arma_abs(val) <= abs_limit) ? eT(0) : val;
144 }
145 }
146
147
148
149 template<typename T>
150 arma_hot
151 inline
152 void
clean(std::complex<T> * mem,const uword n_elem,const T abs_limit)153 arrayops::clean(std::complex<T>* mem, const uword n_elem, const T abs_limit)
154 {
155 typedef typename std::complex<T> eT;
156
157 for(uword i=0; i<n_elem; ++i)
158 {
159 eT& val = mem[i];
160
161 T val_real = std::real(val);
162 T val_imag = std::imag(val);
163
164 if(std::abs(val_real) <= abs_limit)
165 {
166 val_imag = (std::abs(val_imag) <= abs_limit) ? T(0) : val_imag;
167
168 val = std::complex<T>(T(0), val_imag);
169 }
170 else
171 if(std::abs(val_imag) <= abs_limit)
172 {
173 val = std::complex<T>(val_real, T(0));
174 }
175 }
176 }
177
178
179
180 template<typename eT>
181 inline
182 void
clamp(eT * mem,const uword n_elem,const eT min_val,const eT max_val,const typename arma_not_cx<eT>::result * junk)183 arrayops::clamp(eT* mem, const uword n_elem, const eT min_val, const eT max_val, const typename arma_not_cx<eT>::result* junk)
184 {
185 arma_ignore(junk);
186
187 for(uword i=0; i<n_elem; ++i)
188 {
189 eT& val = mem[i];
190
191 val = (val < min_val) ? min_val : ((val > max_val) ? max_val : val);
192 }
193 }
194
195
196
197 template<typename T>
198 inline
199 void
clamp(std::complex<T> * mem,const uword n_elem,const std::complex<T> & min_val,const std::complex<T> & max_val)200 arrayops::clamp(std::complex<T>* mem, const uword n_elem, const std::complex<T>& min_val, const std::complex<T>& max_val)
201 {
202 typedef typename std::complex<T> eT;
203
204 const T min_val_real = std::real(min_val);
205 const T min_val_imag = std::imag(min_val);
206
207 const T max_val_real = std::real(max_val);
208 const T max_val_imag = std::imag(max_val);
209
210 for(uword i=0; i<n_elem; ++i)
211 {
212 eT& val = mem[i];
213
214 T val_real = std::real(val);
215 T val_imag = std::imag(val);
216
217 val_real = (val_real < min_val_real) ? min_val_real : ((val_real > max_val_real) ? max_val_real : val_real);
218 val_imag = (val_imag < min_val_imag) ? min_val_imag : ((val_imag > max_val_imag) ? max_val_imag : val_imag);
219
220 val = std::complex<T>(val_real,val_imag);
221 }
222 }
223
224
225
226 template<typename out_eT, typename in_eT>
227 arma_inline
228 void
convert_cx_scalar(out_eT & out,const in_eT & in,const typename arma_not_cx<out_eT>::result * junk1,const typename arma_not_cx<in_eT>::result * junk2)229 arrayops::convert_cx_scalar
230 (
231 out_eT& out,
232 const in_eT& in,
233 const typename arma_not_cx<out_eT>::result* junk1,
234 const typename arma_not_cx< in_eT>::result* junk2
235 )
236 {
237 arma_ignore(junk1);
238 arma_ignore(junk2);
239
240 out = out_eT(in);
241 }
242
243
244
245 template<typename out_eT, typename in_T>
246 arma_inline
247 void
convert_cx_scalar(out_eT & out,const std::complex<in_T> & in,const typename arma_not_cx<out_eT>::result * junk)248 arrayops::convert_cx_scalar
249 (
250 out_eT& out,
251 const std::complex<in_T>& in,
252 const typename arma_not_cx<out_eT>::result* junk
253 )
254 {
255 arma_ignore(junk);
256
257 const in_T val = in.real();
258
259 const bool conversion_ok = (std::is_integral<out_eT>::value && std::is_floating_point<in_T>::value) ? arma_isfinite(val) : true;
260
261 out = conversion_ok ? out_eT(val) : out_eT(0);
262 }
263
264
265
266 template<typename out_T, typename in_T>
267 arma_inline
268 void
convert_cx_scalar(std::complex<out_T> & out,const std::complex<in_T> & in)269 arrayops::convert_cx_scalar
270 (
271 std::complex<out_T>& out,
272 const std::complex< in_T>& in
273 )
274 {
275 typedef std::complex<out_T> out_eT;
276
277 out = out_eT(in);
278 }
279
280
281
282 template<typename out_eT, typename in_eT>
283 arma_hot
284 inline
285 void
convert(out_eT * dest,const in_eT * src,const uword n_elem)286 arrayops::convert(out_eT* dest, const in_eT* src, const uword n_elem)
287 {
288 if(is_same_type<out_eT,in_eT>::value)
289 {
290 const out_eT* src2 = (const out_eT*)src;
291
292 if(dest != src2) { arrayops::copy(dest, src2, n_elem); }
293
294 return;
295 }
296
297 const bool check_finite = (std::is_integral<out_eT>::value && std::is_floating_point<in_eT>::value);
298
299 uword j;
300
301 for(j=1; j<n_elem; j+=2)
302 {
303 const in_eT tmp_i = (*src); src++;
304 const in_eT tmp_j = (*src); src++;
305
306 // dest[i] = out_eT( tmp_i );
307 // dest[j] = out_eT( tmp_j );
308
309 const bool ok_i = check_finite ? arma_isfinite(tmp_i) : true;
310 const bool ok_j = check_finite ? arma_isfinite(tmp_j) : true;
311
312 (*dest) = ok_i
313 ? (
314 (is_signed<out_eT>::value)
315 ? out_eT( tmp_i )
316 : ( cond_rel< is_signed<in_eT>::value >::lt(tmp_i, in_eT(0)) ? out_eT(0) : out_eT(tmp_i) )
317 )
318 : out_eT(0);
319
320 dest++;
321
322 (*dest) = ok_j
323 ? (
324 (is_signed<out_eT>::value)
325 ? out_eT( tmp_j )
326 : ( cond_rel< is_signed<in_eT>::value >::lt(tmp_j, in_eT(0)) ? out_eT(0) : out_eT(tmp_j) )
327 )
328 : out_eT(0);
329 dest++;
330 }
331
332 if((j-1) < n_elem)
333 {
334 const in_eT tmp_i = (*src);
335
336 // dest[i] = out_eT( tmp_i );
337
338 const bool ok_i = check_finite ? arma_isfinite(tmp_i) : true;
339
340 (*dest) = ok_i
341 ? (
342 (is_signed<out_eT>::value)
343 ? out_eT( tmp_i )
344 : ( cond_rel< is_signed<in_eT>::value >::lt(tmp_i, in_eT(0)) ? out_eT(0) : out_eT(tmp_i) )
345 )
346 : out_eT(0);
347 }
348 }
349
350
351
352 template<typename out_eT, typename in_eT>
353 arma_hot
354 inline
355 void
convert_cx(out_eT * dest,const in_eT * src,const uword n_elem)356 arrayops::convert_cx(out_eT* dest, const in_eT* src, const uword n_elem)
357 {
358 uword j;
359
360 for(j=1; j<n_elem; j+=2)
361 {
362 arrayops::convert_cx_scalar( (*dest), (*src) ); dest++; src++;
363 arrayops::convert_cx_scalar( (*dest), (*src) ); dest++; src++;
364 }
365
366 if((j-1) < n_elem)
367 {
368 arrayops::convert_cx_scalar( (*dest), (*src) );
369 }
370 }
371
372
373
374 template<typename eT>
375 arma_hot
376 inline
377 void
inplace_plus(eT * dest,const eT * src,const uword n_elem)378 arrayops::inplace_plus(eT* dest, const eT* src, const uword n_elem)
379 {
380 if(memory::is_aligned(dest))
381 {
382 memory::mark_as_aligned(dest);
383
384 if(memory::is_aligned(src))
385 {
386 memory::mark_as_aligned(src);
387
388 arrayops::inplace_plus_base(dest, src, n_elem);
389 }
390 else
391 {
392 arrayops::inplace_plus_base(dest, src, n_elem);
393 }
394 }
395 else
396 {
397 if(memory::is_aligned(src))
398 {
399 memory::mark_as_aligned(src);
400
401 arrayops::inplace_plus_base(dest, src, n_elem);
402 }
403 else
404 {
405 arrayops::inplace_plus_base(dest, src, n_elem);
406 }
407 }
408 }
409
410
411
412 template<typename eT>
413 arma_hot
414 inline
415 void
inplace_minus(eT * dest,const eT * src,const uword n_elem)416 arrayops::inplace_minus(eT* dest, const eT* src, const uword n_elem)
417 {
418 if(memory::is_aligned(dest))
419 {
420 memory::mark_as_aligned(dest);
421
422 if(memory::is_aligned(src))
423 {
424 memory::mark_as_aligned(src);
425
426 arrayops::inplace_minus_base(dest, src, n_elem);
427 }
428 else
429 {
430 arrayops::inplace_minus_base(dest, src, n_elem);
431 }
432 }
433 else
434 {
435 if(memory::is_aligned(src))
436 {
437 memory::mark_as_aligned(src);
438
439 arrayops::inplace_minus_base(dest, src, n_elem);
440 }
441 else
442 {
443 arrayops::inplace_minus_base(dest, src, n_elem);
444 }
445 }
446 }
447
448
449
450 template<typename eT>
451 arma_hot
452 inline
453 void
inplace_mul(eT * dest,const eT * src,const uword n_elem)454 arrayops::inplace_mul(eT* dest, const eT* src, const uword n_elem)
455 {
456 if(memory::is_aligned(dest))
457 {
458 memory::mark_as_aligned(dest);
459
460 if(memory::is_aligned(src))
461 {
462 memory::mark_as_aligned(src);
463
464 arrayops::inplace_mul_base(dest, src, n_elem);
465 }
466 else
467 {
468 arrayops::inplace_mul_base(dest, src, n_elem);
469 }
470 }
471 else
472 {
473 if(memory::is_aligned(src))
474 {
475 memory::mark_as_aligned(src);
476
477 arrayops::inplace_mul_base(dest, src, n_elem);
478 }
479 else
480 {
481 arrayops::inplace_mul_base(dest, src, n_elem);
482 }
483 }
484 }
485
486
487
488 template<typename eT>
489 arma_hot
490 inline
491 void
inplace_div(eT * dest,const eT * src,const uword n_elem)492 arrayops::inplace_div(eT* dest, const eT* src, const uword n_elem)
493 {
494 if(memory::is_aligned(dest))
495 {
496 memory::mark_as_aligned(dest);
497
498 if(memory::is_aligned(src))
499 {
500 memory::mark_as_aligned(src);
501
502 arrayops::inplace_div_base(dest, src, n_elem);
503 }
504 else
505 {
506 arrayops::inplace_div_base(dest, src, n_elem);
507 }
508 }
509 else
510 {
511 if(memory::is_aligned(src))
512 {
513 memory::mark_as_aligned(src);
514
515 arrayops::inplace_div_base(dest, src, n_elem);
516 }
517 else
518 {
519 arrayops::inplace_div_base(dest, src, n_elem);
520 }
521 }
522 }
523
524
525
526 template<typename eT>
527 arma_hot
528 inline
529 void
inplace_plus_base(eT * dest,const eT * src,const uword n_elem)530 arrayops::inplace_plus_base(eT* dest, const eT* src, const uword n_elem)
531 {
532 #if defined(ARMA_SIMPLE_LOOPS)
533 {
534 for(uword i=0; i<n_elem; ++i)
535 {
536 dest[i] += src[i];
537 }
538 }
539 #else
540 {
541 uword i,j;
542
543 for(i=0, j=1; j<n_elem; i+=2, j+=2)
544 {
545 const eT tmp_i = src[i];
546 const eT tmp_j = src[j];
547
548 dest[i] += tmp_i;
549 dest[j] += tmp_j;
550 }
551
552 if(i < n_elem)
553 {
554 dest[i] += src[i];
555 }
556 }
557 #endif
558 }
559
560
561
562 template<typename eT>
563 arma_hot
564 inline
565 void
inplace_minus_base(eT * dest,const eT * src,const uword n_elem)566 arrayops::inplace_minus_base(eT* dest, const eT* src, const uword n_elem)
567 {
568 #if defined(ARMA_SIMPLE_LOOPS)
569 {
570 for(uword i=0; i<n_elem; ++i)
571 {
572 dest[i] -= src[i];
573 }
574 }
575 #else
576 {
577 uword i,j;
578
579 for(i=0, j=1; j<n_elem; i+=2, j+=2)
580 {
581 const eT tmp_i = src[i];
582 const eT tmp_j = src[j];
583
584 dest[i] -= tmp_i;
585 dest[j] -= tmp_j;
586 }
587
588 if(i < n_elem)
589 {
590 dest[i] -= src[i];
591 }
592 }
593 #endif
594 }
595
596
597
598 template<typename eT>
599 arma_hot
600 inline
601 void
inplace_mul_base(eT * dest,const eT * src,const uword n_elem)602 arrayops::inplace_mul_base(eT* dest, const eT* src, const uword n_elem)
603 {
604 #if defined(ARMA_SIMPLE_LOOPS)
605 {
606 for(uword i=0; i<n_elem; ++i)
607 {
608 dest[i] *= src[i];
609 }
610 }
611 #else
612 {
613 uword i,j;
614
615 for(i=0, j=1; j<n_elem; i+=2, j+=2)
616 {
617 const eT tmp_i = src[i];
618 const eT tmp_j = src[j];
619
620 dest[i] *= tmp_i;
621 dest[j] *= tmp_j;
622 }
623
624 if(i < n_elem)
625 {
626 dest[i] *= src[i];
627 }
628 }
629 #endif
630 }
631
632
633
634 template<typename eT>
635 arma_hot
636 inline
637 void
inplace_div_base(eT * dest,const eT * src,const uword n_elem)638 arrayops::inplace_div_base(eT* dest, const eT* src, const uword n_elem)
639 {
640 #if defined(ARMA_SIMPLE_LOOPS)
641 {
642 for(uword i=0; i<n_elem; ++i)
643 {
644 dest[i] /= src[i];
645 }
646 }
647 #else
648 {
649 uword i,j;
650
651 for(i=0, j=1; j<n_elem; i+=2, j+=2)
652 {
653 const eT tmp_i = src[i];
654 const eT tmp_j = src[j];
655
656 dest[i] /= tmp_i;
657 dest[j] /= tmp_j;
658 }
659
660 if(i < n_elem)
661 {
662 dest[i] /= src[i];
663 }
664 }
665 #endif
666 }
667
668
669
670 template<typename eT>
671 arma_hot
672 inline
673 void
inplace_set(eT * dest,const eT val,const uword n_elem)674 arrayops::inplace_set(eT* dest, const eT val, const uword n_elem)
675 {
676 if(val == eT(0))
677 {
678 arrayops::fill_zeros(dest, n_elem);
679 }
680 else
681 {
682 if( (n_elem <= 9) && (is_cx<eT>::no) )
683 {
684 arrayops::inplace_set_small(dest, val, n_elem);
685 }
686 else
687 {
688 arrayops::inplace_set_simple(dest, val, n_elem);
689 }
690 }
691 }
692
693
694
695 template<typename eT>
696 arma_hot
697 inline
698 void
inplace_set_simple(eT * dest,const eT val,const uword n_elem)699 arrayops::inplace_set_simple(eT* dest, const eT val, const uword n_elem)
700 {
701 if(memory::is_aligned(dest))
702 {
703 memory::mark_as_aligned(dest);
704
705 arrayops::inplace_set_base(dest, val, n_elem);
706 }
707 else
708 {
709 arrayops::inplace_set_base(dest, val, n_elem);
710 }
711 }
712
713
714
715 template<typename eT>
716 arma_hot
717 inline
718 void
inplace_set_base(eT * dest,const eT val,const uword n_elem)719 arrayops::inplace_set_base(eT* dest, const eT val, const uword n_elem)
720 {
721 #if defined(ARMA_SIMPLE_LOOPS)
722 {
723 for(uword i=0; i<n_elem; ++i)
724 {
725 dest[i] = val;
726 }
727 }
728 #else
729 {
730 uword i,j;
731
732 for(i=0, j=1; j<n_elem; i+=2, j+=2)
733 {
734 dest[i] = val;
735 dest[j] = val;
736 }
737
738 if(i < n_elem)
739 {
740 dest[i] = val;
741 }
742 }
743 #endif
744 }
745
746
747
748 template<typename eT>
749 arma_cold
750 inline
751 void
inplace_set_small(eT * dest,const eT val,const uword n_elem)752 arrayops::inplace_set_small(eT* dest, const eT val, const uword n_elem)
753 {
754 switch(n_elem)
755 {
756 case 9: dest[ 8] = val;
757 // fallthrough
758 case 8: dest[ 7] = val;
759 // fallthrough
760 case 7: dest[ 6] = val;
761 // fallthrough
762 case 6: dest[ 5] = val;
763 // fallthrough
764 case 5: dest[ 4] = val;
765 // fallthrough
766 case 4: dest[ 3] = val;
767 // fallthrough
768 case 3: dest[ 2] = val;
769 // fallthrough
770 case 2: dest[ 1] = val;
771 // fallthrough
772 case 1: dest[ 0] = val;
773 // fallthrough
774 default:;
775 }
776 }
777
778
779
780 template<typename eT, const uword n_elem>
781 arma_hot
782 inline
783 void
inplace_set_fixed(eT * dest,const eT val)784 arrayops::inplace_set_fixed(eT* dest, const eT val)
785 {
786 for(uword i=0; i<n_elem; ++i)
787 {
788 dest[i] = val;
789 }
790 }
791
792
793
794 template<typename eT>
795 arma_hot
796 inline
797 void
inplace_plus(eT * dest,const eT val,const uword n_elem)798 arrayops::inplace_plus(eT* dest, const eT val, const uword n_elem)
799 {
800 if(memory::is_aligned(dest))
801 {
802 memory::mark_as_aligned(dest);
803
804 arrayops::inplace_plus_base(dest, val, n_elem);
805 }
806 else
807 {
808 arrayops::inplace_plus_base(dest, val, n_elem);
809 }
810 }
811
812
813
814 template<typename eT>
815 arma_hot
816 inline
817 void
inplace_minus(eT * dest,const eT val,const uword n_elem)818 arrayops::inplace_minus(eT* dest, const eT val, const uword n_elem)
819 {
820 if(memory::is_aligned(dest))
821 {
822 memory::mark_as_aligned(dest);
823
824 arrayops::inplace_minus_base(dest, val, n_elem);
825 }
826 else
827 {
828 arrayops::inplace_minus_base(dest, val, n_elem);
829 }
830 }
831
832
833
834 template<typename eT>
835 arma_hot
836 inline
837 void
inplace_mul(eT * dest,const eT val,const uword n_elem)838 arrayops::inplace_mul(eT* dest, const eT val, const uword n_elem)
839 {
840 if(memory::is_aligned(dest))
841 {
842 memory::mark_as_aligned(dest);
843
844 arrayops::inplace_mul_base(dest, val, n_elem);
845 }
846 else
847 {
848 arrayops::inplace_mul_base(dest, val, n_elem);
849 }
850 }
851
852
853
854 template<typename eT>
855 arma_hot
856 inline
857 void
inplace_div(eT * dest,const eT val,const uword n_elem)858 arrayops::inplace_div(eT* dest, const eT val, const uword n_elem)
859 {
860 if(memory::is_aligned(dest))
861 {
862 memory::mark_as_aligned(dest);
863
864 arrayops::inplace_div_base(dest, val, n_elem);
865 }
866 else
867 {
868 arrayops::inplace_div_base(dest, val, n_elem);
869 }
870 }
871
872
873
874 template<typename eT>
875 arma_hot
876 inline
877 void
inplace_plus_base(eT * dest,const eT val,const uword n_elem)878 arrayops::inplace_plus_base(eT* dest, const eT val, const uword n_elem)
879 {
880 #if defined(ARMA_SIMPLE_LOOPS)
881 {
882 for(uword i=0; i<n_elem; ++i)
883 {
884 dest[i] += val;
885 }
886 }
887 #else
888 {
889 uword i,j;
890
891 for(i=0, j=1; j<n_elem; i+=2, j+=2)
892 {
893 dest[i] += val;
894 dest[j] += val;
895 }
896
897 if(i < n_elem)
898 {
899 dest[i] += val;
900 }
901 }
902 #endif
903 }
904
905
906
907 template<typename eT>
908 arma_hot
909 inline
910 void
inplace_minus_base(eT * dest,const eT val,const uword n_elem)911 arrayops::inplace_minus_base(eT* dest, const eT val, const uword n_elem)
912 {
913 #if defined(ARMA_SIMPLE_LOOPS)
914 {
915 for(uword i=0; i<n_elem; ++i)
916 {
917 dest[i] -= val;
918 }
919 }
920 #else
921 {
922 uword i,j;
923
924 for(i=0, j=1; j<n_elem; i+=2, j+=2)
925 {
926 dest[i] -= val;
927 dest[j] -= val;
928 }
929
930 if(i < n_elem)
931 {
932 dest[i] -= val;
933 }
934 }
935 #endif
936 }
937
938
939
940 template<typename eT>
941 arma_hot
942 inline
943 void
inplace_mul_base(eT * dest,const eT val,const uword n_elem)944 arrayops::inplace_mul_base(eT* dest, const eT val, const uword n_elem)
945 {
946 #if defined(ARMA_SIMPLE_LOOPS)
947 {
948 for(uword i=0; i<n_elem; ++i)
949 {
950 dest[i] *= val;
951 }
952 }
953 #else
954 {
955 uword i,j;
956
957 for(i=0, j=1; j<n_elem; i+=2, j+=2)
958 {
959 dest[i] *= val;
960 dest[j] *= val;
961 }
962
963 if(i < n_elem)
964 {
965 dest[i] *= val;
966 }
967 }
968 #endif
969 }
970
971
972
973 template<typename eT>
974 arma_hot
975 inline
976 void
inplace_div_base(eT * dest,const eT val,const uword n_elem)977 arrayops::inplace_div_base(eT* dest, const eT val, const uword n_elem)
978 {
979 #if defined(ARMA_SIMPLE_LOOPS)
980 {
981 for(uword i=0; i<n_elem; ++i)
982 {
983 dest[i] /= val;
984 }
985 }
986 #else
987 {
988 uword i,j;
989
990 for(i=0, j=1; j<n_elem; i+=2, j+=2)
991 {
992 dest[i] /= val;
993 dest[j] /= val;
994 }
995
996 if(i < n_elem)
997 {
998 dest[i] /= val;
999 }
1000 }
1001 #endif
1002 }
1003
1004
1005
1006 template<typename eT>
1007 arma_hot
1008 inline
1009 eT
accumulate(const eT * src,const uword n_elem)1010 arrayops::accumulate(const eT* src, const uword n_elem)
1011 {
1012 #if defined(__FINITE_MATH_ONLY__) && (__FINITE_MATH_ONLY__ > 0)
1013 {
1014 eT acc = eT(0);
1015
1016 if(memory::is_aligned(src))
1017 {
1018 memory::mark_as_aligned(src);
1019 for(uword i=0; i<n_elem; ++i) { acc += src[i]; }
1020 }
1021 else
1022 {
1023 for(uword i=0; i<n_elem; ++i) { acc += src[i]; }
1024 }
1025
1026 return acc;
1027 }
1028 #else
1029 {
1030 eT acc1 = eT(0);
1031 eT acc2 = eT(0);
1032
1033 uword j;
1034
1035 for(j=1; j<n_elem; j+=2)
1036 {
1037 acc1 += (*src); src++;
1038 acc2 += (*src); src++;
1039 }
1040
1041 if((j-1) < n_elem)
1042 {
1043 acc1 += (*src);
1044 }
1045
1046 return acc1 + acc2;
1047 }
1048 #endif
1049 }
1050
1051
1052
1053 template<typename eT>
1054 arma_hot
1055 inline
1056 eT
product(const eT * src,const uword n_elem)1057 arrayops::product(const eT* src, const uword n_elem)
1058 {
1059 eT val1 = eT(1);
1060 eT val2 = eT(1);
1061
1062 uword i,j;
1063
1064 for(i=0, j=1; j<n_elem; i+=2, j+=2)
1065 {
1066 val1 *= src[i];
1067 val2 *= src[j];
1068 }
1069
1070 if(i < n_elem)
1071 {
1072 val1 *= src[i];
1073 }
1074
1075 return val1 * val2;
1076 }
1077
1078
1079
1080 template<typename eT>
1081 arma_hot
1082 inline
1083 bool
is_zero(const eT * mem,const uword n_elem,const eT abs_limit,const typename arma_not_cx<eT>::result * junk)1084 arrayops::is_zero(const eT* mem, const uword n_elem, const eT abs_limit, const typename arma_not_cx<eT>::result* junk)
1085 {
1086 arma_ignore(junk);
1087
1088 if(n_elem == 0) { return false; }
1089
1090 if(abs_limit == eT(0))
1091 {
1092 for(uword i=0; i<n_elem; ++i)
1093 {
1094 if(mem[i] != eT(0)) { return false; }
1095 }
1096 }
1097 else
1098 {
1099 for(uword i=0; i<n_elem; ++i)
1100 {
1101 if(eop_aux::arma_abs(mem[i]) > abs_limit) { return false; }
1102 }
1103 }
1104
1105 return true;
1106 }
1107
1108
1109
1110 template<typename T>
1111 arma_hot
1112 inline
1113 bool
is_zero(const std::complex<T> * mem,const uword n_elem,const T abs_limit)1114 arrayops::is_zero(const std::complex<T>* mem, const uword n_elem, const T abs_limit)
1115 {
1116 typedef typename std::complex<T> eT;
1117
1118 if(n_elem == 0) { return false; }
1119
1120 if(abs_limit == T(0))
1121 {
1122 for(uword i=0; i<n_elem; ++i)
1123 {
1124 const eT& val = mem[i];
1125
1126 if(std::real(val) != T(0)) { return false; }
1127 if(std::imag(val) != T(0)) { return false; }
1128 }
1129 }
1130 else
1131 {
1132 for(uword i=0; i<n_elem; ++i)
1133 {
1134 const eT& val = mem[i];
1135
1136 if(std::abs(std::real(val)) > abs_limit) { return false; }
1137 if(std::abs(std::imag(val)) > abs_limit) { return false; }
1138 }
1139 }
1140
1141 return true;
1142 }
1143
1144
1145
1146 template<typename eT>
1147 arma_hot
1148 inline
1149 bool
is_finite(const eT * src,const uword n_elem)1150 arrayops::is_finite(const eT* src, const uword n_elem)
1151 {
1152 uword j;
1153
1154 for(j=1; j<n_elem; j+=2)
1155 {
1156 const eT val_i = (*src); src++;
1157 const eT val_j = (*src); src++;
1158
1159 if( (arma_isfinite(val_i) == false) || (arma_isfinite(val_j) == false) )
1160 {
1161 return false;
1162 }
1163 }
1164
1165 if((j-1) < n_elem)
1166 {
1167 if(arma_isfinite(*src) == false)
1168 {
1169 return false;
1170 }
1171 }
1172
1173 return true;
1174 }
1175
1176
1177
1178 template<typename eT>
1179 arma_hot
1180 inline
1181 bool
has_inf(const eT * src,const uword n_elem)1182 arrayops::has_inf(const eT* src, const uword n_elem)
1183 {
1184 uword j;
1185
1186 for(j=1; j<n_elem; j+=2)
1187 {
1188 const eT val_i = (*src); src++;
1189 const eT val_j = (*src); src++;
1190
1191 if( arma_isinf(val_i) || arma_isinf(val_j) ) { return true; }
1192 }
1193
1194 if((j-1) < n_elem)
1195 {
1196 if(arma_isinf(*src)) { return true; }
1197 }
1198
1199 return false;
1200 }
1201
1202
1203
1204 template<typename eT>
1205 arma_hot
1206 inline
1207 bool
has_nan(const eT * src,const uword n_elem)1208 arrayops::has_nan(const eT* src, const uword n_elem)
1209 {
1210 uword j;
1211
1212 for(j=1; j<n_elem; j+=2)
1213 {
1214 const eT val_i = (*src); src++;
1215 const eT val_j = (*src); src++;
1216
1217 if( arma_isnan(val_i) || arma_isnan(val_j) ) { return true; }
1218 }
1219
1220 if((j-1) < n_elem)
1221 {
1222 if(arma_isnan(*src)) { return true; }
1223 }
1224
1225 return false;
1226 }
1227
1228
1229
1230 //! @}
1231