MOAB
4.9.3pre
|
00001 // This file is part of Eigen, a lightweight C++ template library 00002 // for linear algebra. 00003 // 00004 // Copyright (C) 2008-2015 Gael Guennebaud <[email protected]> 00005 // 00006 // This Source Code Form is subject to the terms of the Mozilla 00007 // Public License v. 2.0. If a copy of the MPL was not distributed 00008 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 00009 00010 #ifndef EIGEN_SPARSEPRODUCT_H 00011 #define EIGEN_SPARSEPRODUCT_H 00012 00013 namespace Eigen { 00014 00026 template<typename Derived> 00027 template<typename OtherDerived> 00028 inline const Product<Derived,OtherDerived,AliasFreeProduct> 00029 SparseMatrixBase<Derived>::operator*(const SparseMatrixBase<OtherDerived> &other) const 00030 { 00031 return Product<Derived,OtherDerived,AliasFreeProduct>(derived(), other.derived()); 00032 } 00033 00034 namespace internal { 00035 00036 // sparse * sparse 00037 template<typename Lhs, typename Rhs, int ProductType> 00038 struct generic_product_impl<Lhs, Rhs, SparseShape, SparseShape, ProductType> 00039 { 00040 template<typename Dest> 00041 static void evalTo(Dest& dst, const Lhs& lhs, const Rhs& rhs) 00042 { 00043 evalTo(dst, lhs, rhs, typename evaluator_traits<Dest>::Shape()); 00044 } 00045 00046 // dense += sparse * sparse 00047 template<typename Dest,typename ActualLhs> 00048 static void addTo(Dest& dst, const ActualLhs& lhs, const Rhs& rhs, int* = typename enable_if<is_same<typename evaluator_traits<Dest>::Shape,DenseShape>::value,int*>::type(0) ) 00049 { 00050 typedef typename nested_eval<ActualLhs,Dynamic>::type LhsNested; 00051 typedef typename nested_eval<Rhs,Dynamic>::type RhsNested; 00052 LhsNested lhsNested(lhs); 00053 RhsNested rhsNested(rhs); 00054 internal::sparse_sparse_to_dense_product_selector<typename remove_all<LhsNested>::type, 00055 typename remove_all<RhsNested>::type, Dest>::run(lhsNested,rhsNested,dst); 00056 } 00057 00058 // dense -= sparse * sparse 00059 template<typename Dest> 00060 static void subTo(Dest& dst, const Lhs& lhs, const Rhs& rhs, int* = typename enable_if<is_same<typename evaluator_traits<Dest>::Shape,DenseShape>::value,int*>::type(0) ) 00061 { 00062 addTo(dst, -lhs, rhs); 00063 } 00064 00065 protected: 00066 00067 // sparse = sparse * sparse 00068 template<typename Dest> 00069 static void evalTo(Dest& dst, const Lhs& lhs, const Rhs& rhs, SparseShape) 00070 { 00071 typedef typename nested_eval<Lhs,Dynamic>::type LhsNested; 00072 typedef typename nested_eval<Rhs,Dynamic>::type RhsNested; 00073 LhsNested lhsNested(lhs); 00074 RhsNested rhsNested(rhs); 00075 internal::conservative_sparse_sparse_product_selector<typename remove_all<LhsNested>::type, 00076 typename remove_all<RhsNested>::type, Dest>::run(lhsNested,rhsNested,dst); 00077 } 00078 00079 // dense = sparse * sparse 00080 template<typename Dest> 00081 static void evalTo(Dest& dst, const Lhs& lhs, const Rhs& rhs, DenseShape) 00082 { 00083 dst.setZero(); 00084 addTo(dst, lhs, rhs); 00085 } 00086 }; 00087 00088 // sparse * sparse-triangular 00089 template<typename Lhs, typename Rhs, int ProductType> 00090 struct generic_product_impl<Lhs, Rhs, SparseShape, SparseTriangularShape, ProductType> 00091 : public generic_product_impl<Lhs, Rhs, SparseShape, SparseShape, ProductType> 00092 {}; 00093 00094 // sparse-triangular * sparse 00095 template<typename Lhs, typename Rhs, int ProductType> 00096 struct generic_product_impl<Lhs, Rhs, SparseTriangularShape, SparseShape, ProductType> 00097 : public generic_product_impl<Lhs, Rhs, SparseShape, SparseShape, ProductType> 00098 {}; 00099 00100 // dense = sparse-product (can be sparse*sparse, sparse*perm, etc.) 00101 template< typename DstXprType, typename Lhs, typename Rhs> 00102 struct Assignment<DstXprType, Product<Lhs,Rhs,AliasFreeProduct>, internal::assign_op<typename DstXprType::Scalar>, Sparse2Dense> 00103 { 00104 typedef Product<Lhs,Rhs,AliasFreeProduct> SrcXprType; 00105 static void run(DstXprType &dst, const SrcXprType &src, const internal::assign_op<typename DstXprType::Scalar> &) 00106 { 00107 generic_product_impl<Lhs, Rhs>::evalTo(dst,src.lhs(),src.rhs()); 00108 } 00109 }; 00110 00111 // dense += sparse-product (can be sparse*sparse, sparse*perm, etc.) 00112 template< typename DstXprType, typename Lhs, typename Rhs> 00113 struct Assignment<DstXprType, Product<Lhs,Rhs,AliasFreeProduct>, internal::add_assign_op<typename DstXprType::Scalar>, Sparse2Dense> 00114 { 00115 typedef Product<Lhs,Rhs,AliasFreeProduct> SrcXprType; 00116 static void run(DstXprType &dst, const SrcXprType &src, const internal::add_assign_op<typename DstXprType::Scalar> &) 00117 { 00118 generic_product_impl<Lhs, Rhs>::addTo(dst,src.lhs(),src.rhs()); 00119 } 00120 }; 00121 00122 // dense -= sparse-product (can be sparse*sparse, sparse*perm, etc.) 00123 template< typename DstXprType, typename Lhs, typename Rhs> 00124 struct Assignment<DstXprType, Product<Lhs,Rhs,AliasFreeProduct>, internal::sub_assign_op<typename DstXprType::Scalar>, Sparse2Dense> 00125 { 00126 typedef Product<Lhs,Rhs,AliasFreeProduct> SrcXprType; 00127 static void run(DstXprType &dst, const SrcXprType &src, const internal::sub_assign_op<typename DstXprType::Scalar> &) 00128 { 00129 generic_product_impl<Lhs, Rhs>::subTo(dst,src.lhs(),src.rhs()); 00130 } 00131 }; 00132 00133 template<typename Lhs, typename Rhs, int Options> 00134 struct evaluator<SparseView<Product<Lhs, Rhs, Options> > > 00135 : public evaluator<typename Product<Lhs, Rhs, DefaultProduct>::PlainObject> 00136 { 00137 typedef SparseView<Product<Lhs, Rhs, Options> > XprType; 00138 typedef typename XprType::PlainObject PlainObject; 00139 typedef evaluator<PlainObject> Base; 00140 00141 explicit evaluator(const XprType& xpr) 00142 : m_result(xpr.rows(), xpr.cols()) 00143 { 00144 using std::abs; 00145 ::new (static_cast<Base*>(this)) Base(m_result); 00146 typedef typename nested_eval<Lhs,Dynamic>::type LhsNested; 00147 typedef typename nested_eval<Rhs,Dynamic>::type RhsNested; 00148 LhsNested lhsNested(xpr.nestedExpression().lhs()); 00149 RhsNested rhsNested(xpr.nestedExpression().rhs()); 00150 00151 internal::sparse_sparse_product_with_pruning_selector<typename remove_all<LhsNested>::type, 00152 typename remove_all<RhsNested>::type, PlainObject>::run(lhsNested,rhsNested,m_result, 00153 abs(xpr.reference())*xpr.epsilon()); 00154 } 00155 00156 protected: 00157 PlainObject m_result; 00158 }; 00159 00160 } // end namespace internal 00161 00162 } // end namespace Eigen 00163 00164 #endif // EIGEN_SPARSEPRODUCT_H