@@ -174,6 +174,31 @@ op_diagmat::apply(Mat<typename T1::elem_type>& out, const Op< Glue<T1,T2,glue_ti
174174 {
175175 arma_debug_sigprint ();
176176
177+ typedef typename T1 ::elem_type eT;
178+
179+ if (X.m .is_alias (out))
180+ {
181+ Mat<eT> tmp;
182+
183+ op_diagmat::apply_times (tmp, X.m .A , X.m .B );
184+
185+ out.steal_mem (tmp);
186+ }
187+ else
188+ {
189+ op_diagmat::apply_times (out, X.m .A , X.m .B );
190+ }
191+ }
192+
193+
194+
195+ template <typename T1 , typename T2 >
196+ inline
197+ void
198+ op_diagmat::apply (Mat_noalias<typename T1 ::elem_type>& out, const Op< Glue<T1 ,T2 ,glue_times>, op_diagmat>& X)
199+ {
200+ arma_debug_sigprint ();
201+
177202 op_diagmat::apply_times (out, X.m .A , X.m .B );
178203 }
179204
@@ -182,7 +207,7 @@ op_diagmat::apply(Mat<typename T1::elem_type>& out, const Op< Glue<T1,T2,glue_ti
182207template <typename T1 , typename T2 >
183208inline
184209void
185- op_diagmat::apply_times (Mat<typename T1 ::elem_type>& actual_out , const T1 & X, const T2 & Y, const typename arma_not_cx<typename T1 ::elem_type>::result* junk)
210+ op_diagmat::apply_times (Mat<typename T1 ::elem_type>& out , const T1 & X, const T2 & Y, const typename arma_not_cx<typename T1 ::elem_type>::result* junk)
186211 {
187212 arma_debug_sigprint ();
188213 arma_ignore (junk);
@@ -218,9 +243,9 @@ op_diagmat::apply_times(Mat<typename T1::elem_type>& actual_out, const T1& X, co
218243 const eT* C_mem = C.memptr ();
219244 const uword N = C.n_elem ;
220245
221- actual_out .zeros (N,N);
246+ out .zeros (N,N);
222247
223- for (uword i=0 ; i<N; ++i) { actual_out .at (i,i) = (use_alpha) ? eT (alpha * C_mem[i]) : eT (C_mem[i]); }
248+ for (uword i=0 ; i<N; ++i) { out .at (i,i) = (use_alpha) ? eT (alpha * C_mem[i]) : eT (C_mem[i]); }
224249
225250 return ;
226251 }
@@ -236,9 +261,9 @@ op_diagmat::apply_times(Mat<typename T1::elem_type>& actual_out, const T1& X, co
236261 const eT* C_mem = C.memptr ();
237262 const uword N = C.n_elem ;
238263
239- actual_out .zeros (N,N);
264+ out .zeros (N,N);
240265
241- for (uword i=0 ; i<N; ++i) { actual_out .at (i,i) = (use_alpha) ? eT (alpha * C_mem[i]) : eT (C_mem[i]); }
266+ for (uword i=0 ; i<N; ++i) { out .at (i,i) = (use_alpha) ? eT (alpha * C_mem[i]) : eT (C_mem[i]); }
242267
243268 return ;
244269 }
@@ -254,9 +279,9 @@ op_diagmat::apply_times(Mat<typename T1::elem_type>& actual_out, const T1& X, co
254279 const eT* C_mem = C.memptr ();
255280 const uword N = C.n_elem ;
256281
257- actual_out .zeros (N,N);
282+ out .zeros (N,N);
258283
259- for (uword i=0 ; i<N; ++i) { actual_out .at (i,i) = (use_alpha) ? eT (alpha * C_mem[i]) : eT (C_mem[i]); }
284+ for (uword i=0 ; i<N; ++i) { out .at (i,i) = (use_alpha) ? eT (alpha * C_mem[i]) : eT (C_mem[i]); }
260285
261286 return ;
262287 }
@@ -272,20 +297,15 @@ op_diagmat::apply_times(Mat<typename T1::elem_type>& actual_out, const T1& X, co
272297 const eT* C_mem = C.memptr ();
273298 const uword N = C.n_elem ;
274299
275- actual_out .zeros (N,N);
300+ out .zeros (N,N);
276301
277- for (uword i=0 ; i<N; ++i) { actual_out .at (i,i) = (use_alpha) ? eT (alpha * C_mem[i]) : eT (C_mem[i]); }
302+ for (uword i=0 ; i<N; ++i) { out .at (i,i) = (use_alpha) ? eT (alpha * C_mem[i]) : eT (C_mem[i]); }
278303
279304 return ;
280305 }
281306 }
282307
283308 // if we got to this point, the multiplication results in a matrix
284-
285- const bool is_alias = (UA .is_alias (actual_out) || UB .is_alias (actual_out));
286-
287- Mat<eT> tmp;
288- Mat<eT>& out = (is_alias) ? tmp : actual_out;
289309
290310 if ( (partial_unwrap<T1 >::do_trans == false ) && (partial_unwrap<T2 >::do_trans == false ) )
291311 {
@@ -398,16 +418,14 @@ op_diagmat::apply_times(Mat<typename T1::elem_type>& actual_out, const T1& X, co
398418 out.at (k,k) = (use_alpha) ? eT (alpha * acc) : eT (acc);
399419 }
400420 }
401-
402- if (is_alias) { actual_out.steal_mem (tmp); }
403421 }
404422
405423
406424
407425template <typename T1 , typename T2 >
408426inline
409427void
410- op_diagmat::apply_times (Mat<typename T1 ::elem_type>& actual_out , const T1 & X, const T2 & Y, const typename arma_cx_only<typename T1 ::elem_type>::result* junk)
428+ op_diagmat::apply_times (Mat<typename T1 ::elem_type>& out , const T1 & X, const T2 & Y, const typename arma_cx_only<typename T1 ::elem_type>::result* junk)
411429 {
412430 arma_debug_sigprint ();
413431 arma_ignore (junk);
@@ -444,9 +462,9 @@ op_diagmat::apply_times(Mat<typename T1::elem_type>& actual_out, const T1& X, co
444462 const eT* C_mem = C.memptr ();
445463 const uword N = C.n_elem ;
446464
447- actual_out .zeros (N,N);
465+ out .zeros (N,N);
448466
449- for (uword i=0 ; i<N; ++i) { actual_out .at (i,i) = (use_alpha) ? eT (alpha * C_mem[i]) : eT (C_mem[i]); }
467+ for (uword i=0 ; i<N; ++i) { out .at (i,i) = (use_alpha) ? eT (alpha * C_mem[i]) : eT (C_mem[i]); }
450468
451469 return ;
452470 }
@@ -462,9 +480,9 @@ op_diagmat::apply_times(Mat<typename T1::elem_type>& actual_out, const T1& X, co
462480 const eT* C_mem = C.memptr ();
463481 const uword N = C.n_elem ;
464482
465- actual_out .zeros (N,N);
483+ out .zeros (N,N);
466484
467- for (uword i=0 ; i<N; ++i) { actual_out .at (i,i) = (use_alpha) ? eT (alpha * C_mem[i]) : eT (C_mem[i]); }
485+ for (uword i=0 ; i<N; ++i) { out .at (i,i) = (use_alpha) ? eT (alpha * C_mem[i]) : eT (C_mem[i]); }
468486
469487 return ;
470488 }
@@ -480,9 +498,9 @@ op_diagmat::apply_times(Mat<typename T1::elem_type>& actual_out, const T1& X, co
480498 const eT* C_mem = C.memptr ();
481499 const uword N = C.n_elem ;
482500
483- actual_out .zeros (N,N);
501+ out .zeros (N,N);
484502
485- for (uword i=0 ; i<N; ++i) { actual_out .at (i,i) = (use_alpha) ? eT (alpha * C_mem[i]) : eT (C_mem[i]); }
503+ for (uword i=0 ; i<N; ++i) { out .at (i,i) = (use_alpha) ? eT (alpha * C_mem[i]) : eT (C_mem[i]); }
486504
487505 return ;
488506 }
@@ -498,20 +516,15 @@ op_diagmat::apply_times(Mat<typename T1::elem_type>& actual_out, const T1& X, co
498516 const eT* C_mem = C.memptr ();
499517 const uword N = C.n_elem ;
500518
501- actual_out .zeros (N,N);
519+ out .zeros (N,N);
502520
503- for (uword i=0 ; i<N; ++i) { actual_out .at (i,i) = (use_alpha) ? eT (alpha * C_mem[i]) : eT (C_mem[i]); }
521+ for (uword i=0 ; i<N; ++i) { out .at (i,i) = (use_alpha) ? eT (alpha * C_mem[i]) : eT (C_mem[i]); }
504522
505523 return ;
506524 }
507525 }
508526
509527 // if we got to this point, the multiplication results in a matrix
510-
511- const bool is_alias = (UA .is_alias (actual_out) || UB .is_alias (actual_out));
512-
513- Mat<eT> tmp;
514- Mat<eT>& out = (is_alias) ? tmp : actual_out;
515528
516529 if ( (partial_unwrap<T1 >::do_trans == false ) && (partial_unwrap<T2 >::do_trans == false ) )
517530 {
@@ -673,8 +686,6 @@ op_diagmat::apply_times(Mat<typename T1::elem_type>& actual_out, const T1& X, co
673686 out.at (k,k) = (use_alpha) ? eT (alpha * acc) : eT (acc);
674687 }
675688 }
676-
677- if (is_alias) { actual_out.steal_mem (tmp); }
678689 }
679690
680691
0 commit comments