Chris@49: // Copyright (C) 2008-2012 NICTA (www.nicta.com.au)
Chris@49: // Copyright (C) 2008-2012 Conrad Sanderson
Chris@49: // Copyright (C) 2012 Ryan Curtin
Chris@49: // 
Chris@49: // This Source Code Form is subject to the terms of the Mozilla Public
Chris@49: // License, v. 2.0. If a copy of the MPL was not distributed with this
Chris@49: // file, You can obtain one at http://mozilla.org/MPL/2.0/.
Chris@49: 
Chris@49: 
Chris@49: //! \addtogroup fn_dot
Chris@49: //! @{
Chris@49: 
Chris@49: 
Chris@49: template<typename T1, typename T2>
Chris@49: arma_inline
Chris@49: arma_warn_unused
Chris@49: typename
Chris@49: enable_if2
Chris@49:   <
Chris@49:   is_arma_type<T1>::value && is_arma_type<T2>::value && is_same_type<typename T1::elem_type, typename T2::elem_type>::value,
Chris@49:   typename T1::elem_type
Chris@49:   >::result
Chris@49: dot
Chris@49:   (
Chris@49:   const T1& A,
Chris@49:   const T2& B
Chris@49:   )
Chris@49:   {
Chris@49:   arma_extra_debug_sigprint();
Chris@49:   
Chris@49:   return op_dot::apply(A,B);
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename T1, typename T2>
Chris@49: arma_inline
Chris@49: arma_warn_unused
Chris@49: typename
Chris@49: enable_if2
Chris@49:   <
Chris@49:   is_arma_type<T1>::value && is_arma_type<T2>::value && is_same_type<typename T1::elem_type, typename T2::elem_type>::value,
Chris@49:   typename T1::elem_type
Chris@49:   >::result
Chris@49: norm_dot
Chris@49:   (
Chris@49:   const T1& A, 
Chris@49:   const T2& B
Chris@49:   )
Chris@49:   {
Chris@49:   arma_extra_debug_sigprint();
Chris@49:   
Chris@49:   return op_norm_dot::apply(A,B);
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: //
Chris@49: // cdot
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename T1, typename T2>
Chris@49: arma_inline
Chris@49: arma_warn_unused
Chris@49: typename
Chris@49: enable_if2
Chris@49:   <
Chris@49:   is_arma_type<T1>::value && is_arma_type<T2>::value && is_same_type<typename T1::elem_type, typename T2::elem_type>::value && is_not_complex<typename T1::elem_type>::value,
Chris@49:   typename T1::elem_type
Chris@49:   >::result
Chris@49: cdot
Chris@49:   (
Chris@49:   const T1& A,
Chris@49:   const T2& B
Chris@49:   )
Chris@49:   {
Chris@49:   arma_extra_debug_sigprint();
Chris@49:   
Chris@49:   return op_dot::apply(A,B);
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: template<typename T1, typename T2>
Chris@49: arma_inline
Chris@49: arma_warn_unused
Chris@49: typename
Chris@49: enable_if2
Chris@49:   <
Chris@49:   is_arma_type<T1>::value && is_arma_type<T2>::value && is_same_type<typename T1::elem_type, typename T2::elem_type>::value && is_complex<typename T1::elem_type>::value,
Chris@49:   typename T1::elem_type
Chris@49:   >::result
Chris@49: cdot
Chris@49:   (
Chris@49:   const T1& A,
Chris@49:   const T2& B
Chris@49:   )
Chris@49:   {
Chris@49:   arma_extra_debug_sigprint();
Chris@49:   
Chris@49:   return op_cdot::apply(A,B);
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: // convert dot(htrans(x), y) to cdot(x,y)
Chris@49: 
Chris@49: template<typename T1, typename T2>
Chris@49: arma_inline
Chris@49: arma_warn_unused
Chris@49: typename
Chris@49: enable_if2
Chris@49:   <
Chris@49:   is_arma_type<T2>::value && is_same_type<typename T1::elem_type, typename T2::elem_type>::value && is_complex<typename T1::elem_type>::value,
Chris@49:   typename T1::elem_type
Chris@49:   >::result
Chris@49: dot
Chris@49:   (
Chris@49:   const Op<T1, op_htrans>& A,
Chris@49:   const T2&                B
Chris@49:   )
Chris@49:   {
Chris@49:   arma_extra_debug_sigprint();
Chris@49:   
Chris@49:   return cdot(A.m, B);
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: //
Chris@49: // for sparse matrices
Chris@49: // 
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: namespace priv
Chris@49:   {
Chris@49:   
Chris@49:   template<typename T1, typename T2>
Chris@49:   arma_hot
Chris@49:   inline
Chris@49:   typename T1::elem_type
Chris@49:   dot_helper(const SpProxy<T1>& pa, const SpProxy<T2>& pb)
Chris@49:     {
Chris@49:     typedef typename T1::elem_type eT;
Chris@49:     
Chris@49:     // Iterate over both objects and see when they are the same
Chris@49:     eT result = eT(0);
Chris@49:     
Chris@49:     typename SpProxy<T1>::const_iterator_type a_it  = pa.begin();
Chris@49:     typename SpProxy<T1>::const_iterator_type a_end = pa.end();
Chris@49:     
Chris@49:     typename SpProxy<T2>::const_iterator_type b_it  = pb.begin();
Chris@49:     typename SpProxy<T2>::const_iterator_type b_end = pb.end();
Chris@49:     
Chris@49:     while((a_it != a_end) && (b_it != b_end))
Chris@49:       {
Chris@49:       if(a_it == b_it)
Chris@49:         {
Chris@49:         result += (*a_it) * (*b_it);
Chris@49:         
Chris@49:         ++a_it;
Chris@49:         ++b_it;
Chris@49:         }
Chris@49:       else if((a_it.col() < b_it.col()) || ((a_it.col() == b_it.col()) && (a_it.row() < b_it.row())))
Chris@49:         {
Chris@49:         // a_it is "behind"
Chris@49:         ++a_it;
Chris@49:         }
Chris@49:       else
Chris@49:         {
Chris@49:         // b_it is "behind"
Chris@49:         ++b_it;
Chris@49:         }
Chris@49:       }
Chris@49:     
Chris@49:     return result;
Chris@49:     }
Chris@49:   
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: //! dot product of two sparse objects
Chris@49: template<typename T1, typename T2>
Chris@49: arma_warn_unused
Chris@49: arma_hot
Chris@49: inline
Chris@49: typename
Chris@49: enable_if2
Chris@49:   <(is_arma_sparse_type<T1>::value) && (is_arma_sparse_type<T2>::value) && (is_same_type<typename T1::elem_type, typename T2::elem_type>::value),
Chris@49:    typename T1::elem_type
Chris@49:   >::result
Chris@49: dot
Chris@49:   (
Chris@49:   const T1& x,
Chris@49:   const T2& y
Chris@49:   )
Chris@49:   {
Chris@49:   arma_extra_debug_sigprint();
Chris@49:   
Chris@49:   const SpProxy<T1> pa(x);
Chris@49:   const SpProxy<T2> pb(y);
Chris@49:   
Chris@49:   arma_debug_assert_same_size(pa.get_n_rows(), pa.get_n_cols(), pb.get_n_rows(), pb.get_n_cols(), "dot()");
Chris@49:   
Chris@49:   typedef typename T1::elem_type eT;
Chris@49:   
Chris@49:   typedef typename SpProxy<T1>::stored_type pa_Q_type;
Chris@49:   typedef typename SpProxy<T2>::stored_type pb_Q_type;
Chris@49:   
Chris@49:   if(
Chris@49:          ( (SpProxy<T1>::must_use_iterator == false) && (SpProxy<T2>::must_use_iterator == false) )
Chris@49:       && ( (is_SpMat<pa_Q_type>::value     == true ) && (is_SpMat<pb_Q_type>::value     == true ) )   
Chris@49:     )
Chris@49:     {
Chris@49:     const unwrap_spmat<pa_Q_type> tmp_a(pa.Q);
Chris@49:     const unwrap_spmat<pb_Q_type> tmp_b(pb.Q);
Chris@49:     
Chris@49:     const SpMat<eT>& A = tmp_a.M;
Chris@49:     const SpMat<eT>& B = tmp_b.M;
Chris@49:     
Chris@49:     if( &A == &B )
Chris@49:       {
Chris@49:       // We can do it directly!
Chris@49:       return op_dot::direct_dot_arma(A.n_nonzero, A.values, A.values);
Chris@49:       }
Chris@49:     else
Chris@49:       {
Chris@49:       return priv::dot_helper(pa,pb);
Chris@49:       }
Chris@49:     }
Chris@49:   else
Chris@49:     {
Chris@49:     return priv::dot_helper(pa,pb);
Chris@49:     }
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: //! dot product of one dense and one sparse object
Chris@49: template<typename T1, typename T2>
Chris@49: arma_warn_unused
Chris@49: arma_hot
Chris@49: inline
Chris@49: typename
Chris@49: enable_if2
Chris@49:   <(is_arma_type<T1>::value) && (is_arma_sparse_type<T2>::value) && (is_same_type<typename T1::elem_type, typename T2::elem_type>::value),
Chris@49:    typename T1::elem_type
Chris@49:   >::result
Chris@49: dot
Chris@49:   (
Chris@49:   const T1& x,
Chris@49:   const T2& y
Chris@49:   )
Chris@49:   {
Chris@49:   arma_extra_debug_sigprint();
Chris@49:   
Chris@49:   const   Proxy<T1> pa(x);
Chris@49:   const SpProxy<T2> pb(y);
Chris@49:   
Chris@49:   arma_debug_assert_same_size(pa.get_n_rows(), pa.get_n_cols(), pb.get_n_rows(), pb.get_n_cols(), "dot()");
Chris@49:   
Chris@49:   typedef typename T1::elem_type eT;
Chris@49:   
Chris@49:   eT result = eT(0);
Chris@49:   
Chris@49:   typename SpProxy<T2>::const_iterator_type it     = pb.begin();
Chris@49:   typename SpProxy<T2>::const_iterator_type it_end = pb.end();
Chris@49:   
Chris@49:   // prefer_at_accessor won't save us operations
Chris@49:   while(it != it_end)
Chris@49:     {
Chris@49:     result += (*it) * pa.at(it.row(), it.col());
Chris@49:     ++it;
Chris@49:     }
Chris@49:   
Chris@49:   return result;
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: //! dot product of one sparse and one dense object
Chris@49: template<typename T1, typename T2>
Chris@49: arma_warn_unused
Chris@49: arma_hot
Chris@49: inline
Chris@49: typename
Chris@49: enable_if2
Chris@49:   <(is_arma_sparse_type<T1>::value) && (is_arma_type<T2>::value) && (is_same_type<typename T1::elem_type, typename T2::elem_type>::value),
Chris@49:    typename T1::elem_type
Chris@49:   >::result
Chris@49: dot
Chris@49:   (
Chris@49:   const T1& x,
Chris@49:   const T2& y
Chris@49:   )
Chris@49:   {
Chris@49:   arma_extra_debug_sigprint();
Chris@49:   
Chris@49:   // this is commutative
Chris@49:   return dot(y, x);
Chris@49:   }
Chris@49: 
Chris@49: 
Chris@49: 
Chris@49: //! @}