blob: 249dabc14588fe33190b039ec034e32461302f2f [file] [log] [blame]
// This file is part of Eigen, a lightweight C++ template library
// for linear algebra.
//
// Copyright (C) 2008-2015 Gael Guennebaud <gael.guennebaud@inria.fr>
//
// This Source Code Form is subject to the terms of the Mozilla
// Public License v. 2.0. If a copy of the MPL was not distributed
// with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
#ifndef EIGEN_SPARSEPRODUCT_H
#define EIGEN_SPARSEPRODUCT_H
// IWYU pragma: private
#include "./InternalHeaderCheck.h"
namespace Eigen {
/** \returns an expression of the product of two sparse matrices.
* By default a conservative product preserving the symbolic non zeros is performed.
* The automatic pruning of the small values can be achieved by calling the pruned() function
* in which case a totally different product algorithm is employed:
* \code
* C = (A*B).pruned(); // suppress numerical zeros (exact)
* C = (A*B).pruned(ref);
* C = (A*B).pruned(ref,epsilon);
* \endcode
* where \c ref is a meaningful non zero reference value.
* */
template <typename Derived>
template <typename OtherDerived>
inline const Product<Derived, OtherDerived, AliasFreeProduct> SparseMatrixBase<Derived>::operator*(
const SparseMatrixBase<OtherDerived>& other) const {
return Product<Derived, OtherDerived, AliasFreeProduct>(derived(), other.derived());
}
namespace internal {
// sparse * sparse
template <typename Lhs, typename Rhs, int ProductType>
struct generic_product_impl<Lhs, Rhs, SparseShape, SparseShape, ProductType> {
template <typename Dest>
static void evalTo(Dest& dst, const Lhs& lhs, const Rhs& rhs) {
evalTo(dst, lhs, rhs, typename evaluator_traits<Dest>::Shape());
}
// dense += sparse * sparse
template <typename Dest, typename ActualLhs>
static void addTo(Dest& dst, const ActualLhs& lhs, const Rhs& rhs,
std::enable_if_t<is_same<typename evaluator_traits<Dest>::Shape, DenseShape>::value, int*>* = 0) {
typedef typename nested_eval<ActualLhs, Dynamic>::type LhsNested;
typedef typename nested_eval<Rhs, Dynamic>::type RhsNested;
LhsNested lhsNested(lhs);
RhsNested rhsNested(rhs);
internal::sparse_sparse_to_dense_product_selector<remove_all_t<LhsNested>, remove_all_t<RhsNested>, Dest>::run(
lhsNested, rhsNested, dst);
}
// dense -= sparse * sparse
template <typename Dest>
static void subTo(Dest& dst, const Lhs& lhs, const Rhs& rhs,
std::enable_if_t<is_same<typename evaluator_traits<Dest>::Shape, DenseShape>::value, int*>* = 0) {
addTo(dst, -lhs, rhs);
}
protected:
// sparse = sparse * sparse
template <typename Dest>
static void evalTo(Dest& dst, const Lhs& lhs, const Rhs& rhs, SparseShape) {
typedef typename nested_eval<Lhs, Dynamic>::type LhsNested;
typedef typename nested_eval<Rhs, Dynamic>::type RhsNested;
LhsNested lhsNested(lhs);
RhsNested rhsNested(rhs);
internal::conservative_sparse_sparse_product_selector<remove_all_t<LhsNested>, remove_all_t<RhsNested>, Dest>::run(
lhsNested, rhsNested, dst);
}
// dense = sparse * sparse
template <typename Dest>
static void evalTo(Dest& dst, const Lhs& lhs, const Rhs& rhs, DenseShape) {
dst.setZero();
addTo(dst, lhs, rhs);
}
};
// sparse * sparse-triangular
template <typename Lhs, typename Rhs, int ProductType>
struct generic_product_impl<Lhs, Rhs, SparseShape, SparseTriangularShape, ProductType>
: public generic_product_impl<Lhs, Rhs, SparseShape, SparseShape, ProductType> {};
// sparse-triangular * sparse
template <typename Lhs, typename Rhs, int ProductType>
struct generic_product_impl<Lhs, Rhs, SparseTriangularShape, SparseShape, ProductType>
: public generic_product_impl<Lhs, Rhs, SparseShape, SparseShape, ProductType> {};
// dense = sparse-product (can be sparse*sparse, sparse*perm, etc.)
template <typename DstXprType, typename Lhs, typename Rhs>
struct Assignment<
DstXprType, Product<Lhs, Rhs, AliasFreeProduct>,
internal::assign_op<typename DstXprType::Scalar, typename Product<Lhs, Rhs, AliasFreeProduct>::Scalar>,
Sparse2Dense> {
typedef Product<Lhs, Rhs, AliasFreeProduct> SrcXprType;
static void run(DstXprType& dst, const SrcXprType& src,
const internal::assign_op<typename DstXprType::Scalar, typename SrcXprType::Scalar>&) {
Index dstRows = src.rows();
Index dstCols = src.cols();
if ((dst.rows() != dstRows) || (dst.cols() != dstCols)) dst.resize(dstRows, dstCols);
generic_product_impl<Lhs, Rhs>::evalTo(dst, src.lhs(), src.rhs());
}
};
// dense += sparse-product (can be sparse*sparse, sparse*perm, etc.)
template <typename DstXprType, typename Lhs, typename Rhs>
struct Assignment<
DstXprType, Product<Lhs, Rhs, AliasFreeProduct>,
internal::add_assign_op<typename DstXprType::Scalar, typename Product<Lhs, Rhs, AliasFreeProduct>::Scalar>,
Sparse2Dense> {
typedef Product<Lhs, Rhs, AliasFreeProduct> SrcXprType;
static void run(DstXprType& dst, const SrcXprType& src,
const internal::add_assign_op<typename DstXprType::Scalar, typename SrcXprType::Scalar>&) {
generic_product_impl<Lhs, Rhs>::addTo(dst, src.lhs(), src.rhs());
}
};
// dense -= sparse-product (can be sparse*sparse, sparse*perm, etc.)
template <typename DstXprType, typename Lhs, typename Rhs>
struct Assignment<
DstXprType, Product<Lhs, Rhs, AliasFreeProduct>,
internal::sub_assign_op<typename DstXprType::Scalar, typename Product<Lhs, Rhs, AliasFreeProduct>::Scalar>,
Sparse2Dense> {
typedef Product<Lhs, Rhs, AliasFreeProduct> SrcXprType;
static void run(DstXprType& dst, const SrcXprType& src,
const internal::sub_assign_op<typename DstXprType::Scalar, typename SrcXprType::Scalar>&) {
generic_product_impl<Lhs, Rhs>::subTo(dst, src.lhs(), src.rhs());
}
};
template <typename Lhs, typename Rhs, int Options>
struct unary_evaluator<SparseView<Product<Lhs, Rhs, Options> >, IteratorBased>
: public evaluator<typename Product<Lhs, Rhs, DefaultProduct>::PlainObject> {
typedef SparseView<Product<Lhs, Rhs, Options> > XprType;
typedef typename XprType::PlainObject PlainObject;
typedef evaluator<PlainObject> Base;
explicit unary_evaluator(const XprType& xpr) : m_result(xpr.rows(), xpr.cols()) {
using std::abs;
internal::construct_at<Base>(this, m_result);
typedef typename nested_eval<Lhs, Dynamic>::type LhsNested;
typedef typename nested_eval<Rhs, Dynamic>::type RhsNested;
LhsNested lhsNested(xpr.nestedExpression().lhs());
RhsNested rhsNested(xpr.nestedExpression().rhs());
internal::sparse_sparse_product_with_pruning_selector<remove_all_t<LhsNested>, remove_all_t<RhsNested>,
PlainObject>::run(lhsNested, rhsNested, m_result,
abs(xpr.reference()) * xpr.epsilon());
}
protected:
PlainObject m_result;
};
} // end namespace internal
// sparse matrix = sparse-product (can be sparse*sparse, sparse*perm, etc.)
template <typename Scalar, int Options_, typename StorageIndex_>
template <typename Lhs, typename Rhs>
SparseMatrix<Scalar, Options_, StorageIndex_>& SparseMatrix<Scalar, Options_, StorageIndex_>::operator=(
const Product<Lhs, Rhs, AliasFreeProduct>& src) {
// std::cout << "in Assignment : " << DstOptions << "\n";
SparseMatrix dst(src.rows(), src.cols());
internal::generic_product_impl<Lhs, Rhs>::evalTo(dst, src.lhs(), src.rhs());
this->swap(dst);
return *this;
}
} // end namespace Eigen
#endif // EIGEN_SPARSEPRODUCT_H