Eigen  3.4.90 (git rev 5a9f66fb35d03a4da9ef8976e67a61b30aa16dcf)
 
Loading...
Searching...
No Matches
SparseAssign.h
1// This file is part of Eigen, a lightweight C++ template library
2// for linear algebra.
3//
4// Copyright (C) 2008-2014 Gael Guennebaud <[email protected]>
5//
6// This Source Code Form is subject to the terms of the Mozilla
7// Public License v. 2.0. If a copy of the MPL was not distributed
8// with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
9
10#ifndef EIGEN_SPARSEASSIGN_H
11#define EIGEN_SPARSEASSIGN_H
12
13// IWYU pragma: private
14#include "./InternalHeaderCheck.h"
15
16namespace Eigen {
17
18template <typename Derived>
19template <typename OtherDerived>
20Derived &SparseMatrixBase<Derived>::operator=(const EigenBase<OtherDerived> &other) {
21 internal::call_assignment_no_alias(derived(), other.derived());
22 return derived();
23}
24
25template <typename Derived>
26template <typename OtherDerived>
27Derived &SparseMatrixBase<Derived>::operator=(const ReturnByValue<OtherDerived> &other) {
28 // TODO use the evaluator mechanism
29 other.evalTo(derived());
30 return derived();
31}
32
33template <typename Derived>
34template <typename OtherDerived>
35inline Derived &SparseMatrixBase<Derived>::operator=(const SparseMatrixBase<OtherDerived> &other) {
36 // by default sparse evaluation do not alias, so we can safely bypass the generic call_assignment routine
37 internal::Assignment<Derived, OtherDerived, internal::assign_op<Scalar, typename OtherDerived::Scalar>>::run(
38 derived(), other.derived(), internal::assign_op<Scalar, typename OtherDerived::Scalar>());
39 return derived();
40}
41
42template <typename Derived>
43inline Derived &SparseMatrixBase<Derived>::operator=(const Derived &other) {
44 internal::call_assignment_no_alias(derived(), other.derived());
45 return derived();
46}
47
48namespace internal {
49
50template <>
51struct storage_kind_to_evaluator_kind<Sparse> {
52 typedef IteratorBased Kind;
53};
54
55template <>
56struct storage_kind_to_shape<Sparse> {
57 typedef SparseShape Shape;
58};
59
60struct Sparse2Sparse {};
61struct Sparse2Dense {};
62
63template <>
64struct AssignmentKind<SparseShape, SparseShape> {
65 typedef Sparse2Sparse Kind;
66};
67template <>
68struct AssignmentKind<SparseShape, SparseTriangularShape> {
69 typedef Sparse2Sparse Kind;
70};
71template <>
72struct AssignmentKind<DenseShape, SparseShape> {
73 typedef Sparse2Dense Kind;
74};
75template <>
76struct AssignmentKind<DenseShape, SparseTriangularShape> {
77 typedef Sparse2Dense Kind;
78};
79
80template <typename DstXprType, typename SrcXprType>
81void assign_sparse_to_sparse(DstXprType &dst, const SrcXprType &src) {
82 typedef typename DstXprType::Scalar Scalar;
83 typedef internal::evaluator<DstXprType> DstEvaluatorType;
84 typedef internal::evaluator<SrcXprType> SrcEvaluatorType;
85
86 SrcEvaluatorType srcEvaluator(src);
87
88 constexpr bool transpose = (DstEvaluatorType::Flags & RowMajorBit) != (SrcEvaluatorType::Flags & RowMajorBit);
89 const Index outerEvaluationSize = (SrcEvaluatorType::Flags & RowMajorBit) ? src.rows() : src.cols();
90
91 Index reserveSize = 0;
92 for (Index j = 0; j < outerEvaluationSize; ++j)
93 for (typename SrcEvaluatorType::InnerIterator it(srcEvaluator, j); it; ++it) reserveSize++;
94
95 if ((!transpose) && src.isRValue()) {
96 // eval without temporary
97 dst.resize(src.rows(), src.cols());
98 dst.setZero();
99 dst.reserve(reserveSize);
100 for (Index j = 0; j < outerEvaluationSize; ++j) {
101 dst.startVec(j);
102 for (typename SrcEvaluatorType::InnerIterator it(srcEvaluator, j); it; ++it) {
103 Scalar v = it.value();
104 dst.insertBackByOuterInner(j, it.index()) = v;
105 }
106 }
107 dst.finalize();
108 } else {
109 // eval through a temporary
110 eigen_assert((((internal::traits<DstXprType>::SupportedAccessPatterns & OuterRandomAccessPattern) ==
111 OuterRandomAccessPattern) ||
112 (!((DstEvaluatorType::Flags & RowMajorBit) != (SrcEvaluatorType::Flags & RowMajorBit)))) &&
113 "the transpose operation is supposed to be handled in SparseMatrix::operator=");
114
115 enum { Flip = (DstEvaluatorType::Flags & RowMajorBit) != (SrcEvaluatorType::Flags & RowMajorBit) };
116
117 DstXprType temp(src.rows(), src.cols());
118
119 temp.reserve(reserveSize);
120 for (Index j = 0; j < outerEvaluationSize; ++j) {
121 temp.startVec(j);
122 for (typename SrcEvaluatorType::InnerIterator it(srcEvaluator, j); it; ++it) {
123 Scalar v = it.value();
124 temp.insertBackByOuterInner(Flip ? it.index() : j, Flip ? j : it.index()) = v;
125 }
126 }
127 temp.finalize();
128
129 dst = temp.markAsRValue();
130 }
131}
132
133// Generic Sparse to Sparse assignment
134template <typename DstXprType, typename SrcXprType, typename Functor>
135struct Assignment<DstXprType, SrcXprType, Functor, Sparse2Sparse> {
136 static void run(DstXprType &dst, const SrcXprType &src,
137 const internal::assign_op<typename DstXprType::Scalar, typename SrcXprType::Scalar> & /*func*/) {
138 assign_sparse_to_sparse(dst.derived(), src.derived());
139 }
140};
141
142// Generic Sparse to Dense assignment
143template <typename DstXprType, typename SrcXprType, typename Functor, typename Weak>
144struct Assignment<DstXprType, SrcXprType, Functor, Sparse2Dense, Weak> {
145 static void run(DstXprType &dst, const SrcXprType &src, const Functor &func) {
146 if (internal::is_same<Functor,
147 internal::assign_op<typename DstXprType::Scalar, typename SrcXprType::Scalar>>::value)
148 dst.setZero();
149
150 internal::evaluator<SrcXprType> srcEval(src);
151 resize_if_allowed(dst, src, func);
152 internal::evaluator<DstXprType> dstEval(dst);
153
154 const Index outerEvaluationSize = (internal::evaluator<SrcXprType>::Flags & RowMajorBit) ? src.rows() : src.cols();
155 for (Index j = 0; j < outerEvaluationSize; ++j)
156 for (typename internal::evaluator<SrcXprType>::InnerIterator i(srcEval, j); i; ++i)
157 func.assignCoeff(dstEval.coeffRef(i.row(), i.col()), i.value());
158 }
159};
160
161// Specialization for dense ?= dense +/- sparse and dense ?= sparse +/- dense
162template <typename DstXprType, typename Func1, typename Func2>
163struct assignment_from_dense_op_sparse {
164 template <typename SrcXprType, typename InitialFunc>
165 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void run(DstXprType &dst, const SrcXprType &src,
166 const InitialFunc & /*func*/) {
167#ifdef EIGEN_SPARSE_ASSIGNMENT_FROM_DENSE_OP_SPARSE_PLUGIN
168 EIGEN_SPARSE_ASSIGNMENT_FROM_DENSE_OP_SPARSE_PLUGIN
169#endif
170
171 call_assignment_no_alias(dst, src.lhs(), Func1());
172 call_assignment_no_alias(dst, src.rhs(), Func2());
173 }
174
175 // Specialization for dense1 = sparse + dense2; -> dense1 = dense2; dense1 += sparse;
176 template <typename Lhs, typename Rhs, typename Scalar>
177 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
178 std::enable_if_t<internal::is_same<typename internal::evaluator_traits<Rhs>::Shape, DenseShape>::value>
179 run(DstXprType &dst, const CwiseBinaryOp<internal::scalar_sum_op<Scalar, Scalar>, const Lhs, const Rhs> &src,
180 const internal::assign_op<typename DstXprType::Scalar, Scalar> & /*func*/) {
181#ifdef EIGEN_SPARSE_ASSIGNMENT_FROM_SPARSE_ADD_DENSE_PLUGIN
182 EIGEN_SPARSE_ASSIGNMENT_FROM_SPARSE_ADD_DENSE_PLUGIN
183#endif
184
185 // Apply the dense matrix first, then the sparse one.
186 call_assignment_no_alias(dst, src.rhs(), Func1());
187 call_assignment_no_alias(dst, src.lhs(), Func2());
188 }
189
190 // Specialization for dense1 = sparse - dense2; -> dense1 = -dense2; dense1 += sparse;
191 template <typename Lhs, typename Rhs, typename Scalar>
192 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
193 std::enable_if_t<internal::is_same<typename internal::evaluator_traits<Rhs>::Shape, DenseShape>::value>
194 run(DstXprType &dst,
195 const CwiseBinaryOp<internal::scalar_difference_op<Scalar, Scalar>, const Lhs, const Rhs> &src,
196 const internal::assign_op<typename DstXprType::Scalar, Scalar> & /*func*/) {
197#ifdef EIGEN_SPARSE_ASSIGNMENT_FROM_SPARSE_SUB_DENSE_PLUGIN
198 EIGEN_SPARSE_ASSIGNMENT_FROM_SPARSE_SUB_DENSE_PLUGIN
199#endif
200
201 // Apply the dense matrix first, then the sparse one.
202 call_assignment_no_alias(dst, -src.rhs(), Func1());
203 call_assignment_no_alias(dst, src.lhs(), add_assign_op<typename DstXprType::Scalar, typename Lhs::Scalar>());
204 }
205};
206
207#define EIGEN_CATCH_ASSIGN_DENSE_OP_SPARSE(ASSIGN_OP, BINOP, ASSIGN_OP2) \
208 template <typename DstXprType, typename Lhs, typename Rhs, typename Scalar> \
209 struct Assignment< \
210 DstXprType, CwiseBinaryOp<internal::BINOP<Scalar, Scalar>, const Lhs, const Rhs>, \
211 internal::ASSIGN_OP<typename DstXprType::Scalar, Scalar>, Sparse2Dense, \
212 std::enable_if_t<internal::is_same<typename internal::evaluator_traits<Lhs>::Shape, DenseShape>::value || \
213 internal::is_same<typename internal::evaluator_traits<Rhs>::Shape, DenseShape>::value>> \
214 : assignment_from_dense_op_sparse<DstXprType, \
215 internal::ASSIGN_OP<typename DstXprType::Scalar, typename Lhs::Scalar>, \
216 internal::ASSIGN_OP2<typename DstXprType::Scalar, typename Rhs::Scalar>> {}
217
218EIGEN_CATCH_ASSIGN_DENSE_OP_SPARSE(assign_op, scalar_sum_op, add_assign_op);
219EIGEN_CATCH_ASSIGN_DENSE_OP_SPARSE(add_assign_op, scalar_sum_op, add_assign_op);
220EIGEN_CATCH_ASSIGN_DENSE_OP_SPARSE(sub_assign_op, scalar_sum_op, sub_assign_op);
221
222EIGEN_CATCH_ASSIGN_DENSE_OP_SPARSE(assign_op, scalar_difference_op, sub_assign_op);
223EIGEN_CATCH_ASSIGN_DENSE_OP_SPARSE(add_assign_op, scalar_difference_op, sub_assign_op);
224EIGEN_CATCH_ASSIGN_DENSE_OP_SPARSE(sub_assign_op, scalar_difference_op, add_assign_op);
225
226// Specialization for "dst = dec.solve(rhs)"
227// NOTE we need to specialize it for Sparse2Sparse to avoid ambiguous specialization error
228template <typename DstXprType, typename DecType, typename RhsType, typename Scalar>
229struct Assignment<DstXprType, Solve<DecType, RhsType>, internal::assign_op<Scalar, Scalar>, Sparse2Sparse> {
230 typedef Solve<DecType, RhsType> SrcXprType;
231 static void run(DstXprType &dst, const SrcXprType &src, const internal::assign_op<Scalar, Scalar> &) {
232 Index dstRows = src.rows();
233 Index dstCols = src.cols();
234 if ((dst.rows() != dstRows) || (dst.cols() != dstCols)) dst.resize(dstRows, dstCols);
235
236 src.dec()._solve_impl(src.rhs(), dst);
237 }
238};
239
240struct Diagonal2Sparse {};
241
242template <>
243struct AssignmentKind<SparseShape, DiagonalShape> {
244 typedef Diagonal2Sparse Kind;
245};
246
247template <typename DstXprType, typename SrcXprType, typename Functor>
248struct Assignment<DstXprType, SrcXprType, Functor, Diagonal2Sparse> {
249 typedef typename DstXprType::StorageIndex StorageIndex;
250 typedef typename DstXprType::Scalar Scalar;
251
252 template <int Options, typename AssignFunc>
253 static void run(SparseMatrix<Scalar, Options, StorageIndex> &dst, const SrcXprType &src, const AssignFunc &func) {
254 dst.assignDiagonal(src.diagonal(), func);
255 }
256
257 template <typename DstDerived>
258 static void run(SparseMatrixBase<DstDerived> &dst, const SrcXprType &src,
259 const internal::assign_op<typename DstXprType::Scalar, typename SrcXprType::Scalar> & /*func*/) {
260 dst.derived().diagonal() = src.diagonal();
261 }
262
263 template <typename DstDerived>
264 static void run(SparseMatrixBase<DstDerived> &dst, const SrcXprType &src,
265 const internal::add_assign_op<typename DstXprType::Scalar, typename SrcXprType::Scalar> & /*func*/) {
266 dst.derived().diagonal() += src.diagonal();
267 }
268
269 template <typename DstDerived>
270 static void run(SparseMatrixBase<DstDerived> &dst, const SrcXprType &src,
271 const internal::sub_assign_op<typename DstXprType::Scalar, typename SrcXprType::Scalar> & /*func*/) {
272 dst.derived().diagonal() -= src.diagonal();
273 }
274};
275} // end namespace internal
276
277} // end namespace Eigen
278
279#endif // EIGEN_SPARSEASSIGN_H
Eigen::RowMajorBit
const unsigned int RowMajorBit
Definition Constants.h:70
Eigen
Namespace containing all symbols from the Eigen library.
Definition Core:137