Add hcat/vcat DenseBase concatenation expressions libeigen/eigen!2331
diff --git a/Eigen/Core b/Eigen/Core index 653c266..424396f 100644 --- a/Eigen/Core +++ b/Eigen/Core
@@ -442,6 +442,7 @@ #include "src/Core/PartialReduxEvaluator.h" #include "src/Core/Random.h" #include "src/Core/Replicate.h" +#include "src/Core/ConcatOp.h" #include "src/Core/Reverse.h" #include "src/Core/ArrayWrapper.h" #include "src/Core/StlIterators.h"
diff --git a/Eigen/src/Core/ConcatOp.h b/Eigen/src/Core/ConcatOp.h new file mode 100644 index 0000000..ecf10fa --- /dev/null +++ b/Eigen/src/Core/ConcatOp.h
@@ -0,0 +1,296 @@ +// This file is part of Eigen, a lightweight C++ template library +// for linear algebra. +// +// Copyright (C) 2026 Pavel Guzenfeld +// +// 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_CONCAT_OP_H +#define EIGEN_CONCAT_OP_H + +// IWYU pragma: private +#include "./InternalHeaderCheck.h" + +namespace Eigen { + +namespace internal { + +template <int Direction, typename LhsType, typename RhsType> +struct traits<Concat<Direction, LhsType, RhsType>> : traits<LhsType> { + typedef typename LhsType::Scalar Scalar; + typedef typename traits<LhsType>::StorageKind StorageKind; + typedef typename traits<LhsType>::XprKind XprKind; + typedef typename ref_selector<LhsType>::type LhsTypeNested; + typedef typename ref_selector<RhsType>::type RhsTypeNested; + typedef std::remove_reference_t<LhsTypeNested> LhsTypeNested_; + typedef std::remove_reference_t<RhsTypeNested> RhsTypeNested_; + enum { + // For vertical concat (stacking rows): rows add up, cols must match + // For horizontal concat (stacking cols): cols add up, rows must match + LhsRows = int(LhsType::RowsAtCompileTime), + RhsRows = int(RhsType::RowsAtCompileTime), + LhsCols = int(LhsType::ColsAtCompileTime), + RhsCols = int(RhsType::ColsAtCompileTime), + + RowsAtCompileTime = Direction == Vertical + ? (LhsRows == Dynamic || RhsRows == Dynamic ? int(Dynamic) : LhsRows + RhsRows) + : size_prefer_fixed(LhsRows, RhsRows), + ColsAtCompileTime = Direction == Horizontal + ? (LhsCols == Dynamic || RhsCols == Dynamic ? int(Dynamic) : LhsCols + RhsCols) + : size_prefer_fixed(LhsCols, RhsCols), + + LhsMaxRows = int(LhsType::MaxRowsAtCompileTime), + RhsMaxRows = int(RhsType::MaxRowsAtCompileTime), + LhsMaxCols = int(LhsType::MaxColsAtCompileTime), + RhsMaxCols = int(RhsType::MaxColsAtCompileTime), + + MaxRowsAtCompileTime = + Direction == Vertical + ? (LhsMaxRows == Dynamic || RhsMaxRows == Dynamic ? int(Dynamic) : LhsMaxRows + RhsMaxRows) + : max_size_prefer_dynamic(LhsMaxRows, RhsMaxRows), + MaxColsAtCompileTime = + Direction == Horizontal + ? (LhsMaxCols == Dynamic || RhsMaxCols == Dynamic ? int(Dynamic) : LhsMaxCols + RhsMaxCols) + : max_size_prefer_dynamic(LhsMaxCols, RhsMaxCols), + + IsRowMajor = MaxRowsAtCompileTime == 1 && MaxColsAtCompileTime != 1 ? 1 + : MaxColsAtCompileTime == 1 && MaxRowsAtCompileTime != 1 ? 0 + : (int(LhsType::Flags) & RowMajorBit) ? 1 + : 0, + Flags = IsRowMajor ? RowMajorBit : 0 + }; +}; + +} // namespace internal + +/** + * \class Concat + * \ingroup Core_Module + * + * \brief Expression of the concatenation of two dense expressions + * + * \tparam Direction either \c Vertical or \c Horizontal + * \tparam LhsType the type of the left-hand side expression + * \tparam RhsType the type of the right-hand side expression + * + * This class represents an expression of the concatenation of two dense expressions, + * either vertically (stacking rows) or horizontally (stacking columns). + * + * It is the return type of hcat() and vcat() and typically this is the only way it is used. + * + * \sa hcat(), vcat() + */ +template <int Direction, typename LhsType, typename RhsType> +class Concat : public internal::dense_xpr_base<Concat<Direction, LhsType, RhsType>>::type { + typedef typename internal::traits<Concat>::LhsTypeNested LhsTypeNested; + typedef typename internal::traits<Concat>::RhsTypeNested RhsTypeNested; + typedef typename internal::traits<Concat>::LhsTypeNested_ LhsTypeNested_; + typedef typename internal::traits<Concat>::RhsTypeNested_ RhsTypeNested_; + + public: + typedef typename internal::dense_xpr_base<Concat>::type Base; + EIGEN_DENSE_PUBLIC_INTERFACE(Concat) + typedef internal::remove_all_t<LhsType> LhsNestedExpression; + typedef internal::remove_all_t<RhsType> RhsNestedExpression; + + template <typename OriginalLhsType, typename OriginalRhsType> + EIGEN_DEVICE_FUNC constexpr inline Concat(const OriginalLhsType& lhs, const OriginalRhsType& rhs) + : m_lhs(lhs), m_rhs(rhs) { + EIGEN_STATIC_ASSERT((internal::is_same<std::remove_const_t<LhsType>, OriginalLhsType>::value), + THE_MATRIX_OR_EXPRESSION_THAT_YOU_PASSED_DOES_NOT_HAVE_THE_EXPECTED_TYPE) + EIGEN_STATIC_ASSERT((internal::is_same<std::remove_const_t<RhsType>, OriginalRhsType>::value), + THE_MATRIX_OR_EXPRESSION_THAT_YOU_PASSED_DOES_NOT_HAVE_THE_EXPECTED_TYPE) + EIGEN_STATIC_ASSERT( + (internal::is_same<typename LhsType::Scalar, typename RhsType::Scalar>::value), + YOU_MIXED_DIFFERENT_NUMERIC_TYPES__YOU_NEED_TO_USE_THE_CAST_METHOD_OF_MATRIXBASE_TO_CAST_NUMERIC_TYPES_EXPLICITLY) + EIGEN_STATIC_ASSERT_SAME_XPR_KIND(LhsType, RhsType) + EIGEN_STATIC_ASSERT(Direction != Horizontal || int(LhsType::RowsAtCompileTime) == Dynamic || + int(RhsType::RowsAtCompileTime) == Dynamic || + int(LhsType::RowsAtCompileTime) == int(RhsType::RowsAtCompileTime), + YOU_MIXED_MATRICES_OF_DIFFERENT_SIZES) + EIGEN_STATIC_ASSERT(Direction != Vertical || int(LhsType::ColsAtCompileTime) == Dynamic || + int(RhsType::ColsAtCompileTime) == Dynamic || + int(LhsType::ColsAtCompileTime) == int(RhsType::ColsAtCompileTime), + YOU_MIXED_MATRICES_OF_DIFFERENT_SIZES) + if (Direction == Vertical) { + eigen_assert(lhs.cols() == rhs.cols() && "vcat: number of columns must match"); + } else { + eigen_assert(lhs.rows() == rhs.rows() && "hcat: number of rows must match"); + } + } + + EIGEN_DEVICE_FUNC constexpr Index rows() const { + return Direction == Vertical ? m_lhs.rows() + m_rhs.rows() : m_lhs.rows(); + } + EIGEN_DEVICE_FUNC constexpr Index cols() const { + return Direction == Horizontal ? m_lhs.cols() + m_rhs.cols() : m_lhs.cols(); + } + + EIGEN_DEVICE_FUNC constexpr const LhsTypeNested_& lhs() const { return m_lhs; } + EIGEN_DEVICE_FUNC constexpr const RhsTypeNested_& rhs() const { return m_rhs; } + + protected: + LhsTypeNested m_lhs; + RhsTypeNested m_rhs; +}; + +// Evaluator for Concat +namespace internal { + +template <int Direction, typename LhsType, typename RhsType> +struct evaluator<Concat<Direction, LhsType, RhsType>> : evaluator_base<Concat<Direction, LhsType, RhsType>> { + typedef Concat<Direction, LhsType, RhsType> XprType; + typedef typename XprType::CoeffReturnType CoeffReturnType; + + typedef typename nested_eval<LhsType, 1>::type LhsNested; + typedef typename nested_eval<RhsType, 1>::type RhsNested; + typedef remove_all_t<LhsNested> LhsNestedCleaned; + typedef remove_all_t<RhsNested> RhsNestedCleaned; + + enum { + CoeffReadCost = plain_enum_max(evaluator<LhsNestedCleaned>::CoeffReadCost, + evaluator<RhsNestedCleaned>::CoeffReadCost) + + NumTraits<typename XprType::Scalar>::AddCost, // cost of the branch + LhsFlags = evaluator<LhsNestedCleaned>::Flags, + RhsFlags = evaluator<RhsNestedCleaned>::Flags, + BothHavePacketAccess = (int(LhsFlags) & PacketAccessBit) && (int(RhsFlags) & PacketAccessBit), + BothHaveLinearAccess = (int(LhsFlags) & LinearAccessBit) && (int(RhsFlags) & LinearAccessBit), + IsRowMajor = int(traits<XprType>::Flags) & RowMajorBit, + IsVectorAtCompileTime = XprType::IsVectorAtCompileTime, + Flags = (traits<XprType>::Flags & RowMajorBit) | (BothHavePacketAccess ? PacketAccessBit : 0) | + (IsVectorAtCompileTime && BothHaveLinearAccess ? LinearAccessBit : 0), + Alignment = 0 // conservative: no alignment guarantees across boundary + }; + + EIGEN_DEVICE_FUNC constexpr EIGEN_STRONG_INLINE explicit evaluator(const XprType& xpr) + : m_lhs(xpr.lhs()), + m_rhs(xpr.rhs()), + m_lhsImpl(m_lhs), + m_rhsImpl(m_rhs), + m_lhsRows(xpr.lhs().rows()), + m_lhsCols(xpr.lhs().cols()) {} + + EIGEN_DEVICE_FUNC constexpr EIGEN_STRONG_INLINE CoeffReturnType coeff(Index row, Index col) const { + if (Direction == Vertical) { + if (row < m_lhsRows.value()) + return m_lhsImpl.coeff(row, col); + else + return m_rhsImpl.coeff(row - m_lhsRows.value(), col); + } else { + if (col < m_lhsCols.value()) + return m_lhsImpl.coeff(row, col); + else + return m_rhsImpl.coeff(row, col - m_lhsCols.value()); + } + } + + EIGEN_DEVICE_FUNC constexpr EIGEN_STRONG_INLINE CoeffReturnType coeff(Index index) const { + const Index boundary = Direction == Vertical ? m_lhsRows.value() : m_lhsCols.value(); + if (index < boundary) + return m_lhsImpl.coeff(index); + else + return m_rhsImpl.coeff(index - boundary); + } + + template <int LoadMode, typename PacketType> + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE PacketType packet(Index row, Index col) const { + constexpr int packetSize = unpacket_traits<PacketType>::size; + if (Direction == Vertical) { + const Index boundary = m_lhsRows.value(); + if (row >= boundary) return m_rhsImpl.template packet<LoadMode, PacketType>(row - boundary, col); + // Column-major: inner=rows, packet extends along rows and may straddle the row boundary. + // Row-major: inner=cols, packet extends along cols — never crosses the row boundary. + if (!IsRowMajor && row + packetSize > boundary) return packetBoundary<LoadMode, PacketType>(row, col); + return m_lhsImpl.template packet<LoadMode, PacketType>(row, col); + } else { + const Index boundary = m_lhsCols.value(); + if (col >= boundary) return m_rhsImpl.template packet<LoadMode, PacketType>(row, col - boundary); + // Row-major: inner=cols, packet extends along cols and may straddle the col boundary. + // Column-major: inner=rows, packet extends along rows — never crosses the col boundary. + if (IsRowMajor && col + packetSize > boundary) return packetBoundary<LoadMode, PacketType>(row, col); + return m_lhsImpl.template packet<LoadMode, PacketType>(row, col); + } + } + + template <int LoadMode, typename PacketType> + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE PacketType packet(Index index) const { + constexpr int packetSize = unpacket_traits<PacketType>::size; + const Index boundary = Direction == Vertical ? m_lhsRows.value() : m_lhsCols.value(); + if (index >= boundary) return m_rhsImpl.template packet<LoadMode, PacketType>(index - boundary); + if (index + packetSize > boundary) return packetBoundaryLinear<LoadMode, PacketType>(index); + return m_lhsImpl.template packet<LoadMode, PacketType>(index); + } + + protected: + typedef typename XprType::Scalar Scalar; + + template <int LoadMode, typename PacketType> + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE PacketType packetBoundary(Index row, Index col) const { + constexpr int packetSize = unpacket_traits<PacketType>::size; + EIGEN_ALIGN_MAX Scalar tmp[packetSize]; + for (int i = 0; i < packetSize; ++i) + tmp[i] = coeff(row + (Direction == Vertical ? i : 0), col + (Direction == Horizontal ? i : 0)); + return pload<PacketType>(tmp); + } + + template <int LoadMode, typename PacketType> + EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE PacketType packetBoundaryLinear(Index index) const { + constexpr int packetSize = unpacket_traits<PacketType>::size; + EIGEN_ALIGN_MAX Scalar tmp[packetSize]; + for (int i = 0; i < packetSize; ++i) tmp[i] = coeff(index + i); + return pload<PacketType>(tmp); + } + + const LhsNested m_lhs; + const RhsNested m_rhs; + evaluator<LhsNestedCleaned> m_lhsImpl; + evaluator<RhsNestedCleaned> m_rhsImpl; + const variable_if_dynamic<Index, LhsType::RowsAtCompileTime> m_lhsRows; + const variable_if_dynamic<Index, LhsType::ColsAtCompileTime> m_lhsCols; +}; + +} // namespace internal + +/** + * \relates Concat + * \returns an expression of \a lhs and \a rhs concatenated horizontally (side by side). + * + * Both arguments must have the same number of rows. + * To concatenate more than two expressions, chain calls: \c hcat(hcat(a, b), c). + * + * Example: \code + * Matrix2d A, B; + * auto C = hcat(A, B); // C is 2x4 + * \endcode + * + * \sa vcat(), Concat + */ +template <typename Lhs, typename Rhs> +EIGEN_DEVICE_FUNC inline const Concat<Horizontal, Lhs, Rhs> hcat(const DenseBase<Lhs>& lhs, const DenseBase<Rhs>& rhs) { + return Concat<Horizontal, Lhs, Rhs>(lhs.derived(), rhs.derived()); +} + +/** + * \relates Concat + * \returns an expression of \a lhs and \a rhs concatenated vertically (stacked on top of each other). + * + * Both arguments must have the same number of columns. + * To concatenate more than two expressions, chain calls: \c vcat(vcat(a, b), c). + * + * Example: \code + * Matrix2d A, B; + * auto C = vcat(A, B); // C is 4x2 + * \endcode + * + * \sa hcat(), Concat + */ +template <typename Lhs, typename Rhs> +EIGEN_DEVICE_FUNC inline const Concat<Vertical, Lhs, Rhs> vcat(const DenseBase<Lhs>& lhs, const DenseBase<Rhs>& rhs) { + return Concat<Vertical, Lhs, Rhs>(lhs.derived(), rhs.derived()); +} + +} // end namespace Eigen + +#endif // EIGEN_CONCAT_OP_H
diff --git a/Eigen/src/Core/util/ForwardDeclarations.h b/Eigen/src/Core/util/ForwardDeclarations.h index e2e7155..233e58f 100644 --- a/Eigen/src/Core/util/ForwardDeclarations.h +++ b/Eigen/src/Core/util/ForwardDeclarations.h
@@ -405,6 +405,8 @@ class VectorwiseOp; template <typename MatrixType, int RowFactor, int ColFactor> class Replicate; +template <int Direction, typename LhsType, typename RhsType> +class Concat; template <typename MatrixType, int Direction = BothDirections> class Reverse;
diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index fad0267..b362f1f 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt
@@ -241,6 +241,7 @@ ei_add_test(matrix_cwise) ei_add_test(array_for_matrix) ei_add_test(array_replicate) +ei_add_test(concat) ei_add_test(array_reverse) ei_add_test(ref) ei_add_test(is_same_dense)
diff --git a/test/concat.cpp b/test/concat.cpp new file mode 100644 index 0000000..2f1268d --- /dev/null +++ b/test/concat.cpp
@@ -0,0 +1,854 @@ +// This file is part of Eigen, a lightweight C++ template library +// for linear algebra. +// +// Copyright (C) 2026 Pavel Guzenfeld +// +// 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/. + +#include "main.h" + +using namespace Eigen; + +// ============================================================================ +// Test 1: Dynamic-size matrix concatenation (basic correctness) +// ============================================================================ +template <typename MatrixType> +void test_concat_dynamic(const MatrixType& m) { + typedef typename MatrixType::Scalar Scalar; + typedef Matrix<Scalar, Dynamic, Dynamic> MatrixX; + + Index rows = m.rows(); + Index cols = m.cols(); + + MatrixType a = MatrixType::Random(rows, cols); + MatrixType b = MatrixType::Random(rows, cols); + + // Vertical concatenation: stack rows + { + MatrixX expected(2 * rows, cols); + expected.topRows(rows) = a; + expected.bottomRows(rows) = b; + VERIFY_IS_APPROX(expected, vcat(a, b).eval()); + + // Also verify through assignment to MatrixX + MatrixX result = vcat(a, b); + VERIFY_IS_APPROX(expected, result); + } + + // Horizontal concatenation: stack columns + { + MatrixX expected(rows, 2 * cols); + expected.leftCols(cols) = a; + expected.rightCols(cols) = b; + VERIFY_IS_APPROX(expected, hcat(a, b).eval()); + + MatrixX result = hcat(a, b); + VERIFY_IS_APPROX(expected, result); + } + + // Test with different-sized operands + { + MatrixX c = MatrixX::Random(rows + 2, cols); + MatrixX vresult = vcat(a, c); + VERIFY_IS_EQUAL(vresult.rows(), 2 * rows + 2); + VERIFY_IS_EQUAL(vresult.cols(), cols); + MatrixX vexpected(2 * rows + 2, cols); + vexpected.topRows(rows) = a; + vexpected.bottomRows(rows + 2) = c; + VERIFY_IS_APPROX(vexpected, vresult); + } + + { + MatrixX d = MatrixX::Random(rows, cols + 3); + MatrixX hresult = hcat(a, d); + VERIFY_IS_EQUAL(hresult.rows(), rows); + VERIFY_IS_EQUAL(hresult.cols(), 2 * cols + 3); + MatrixX hexpected(rows, 2 * cols + 3); + hexpected.leftCols(cols) = a; + hexpected.rightCols(cols + 3) = d; + VERIFY_IS_APPROX(hexpected, hresult); + } +} + +// ============================================================================ +// Test 2: Expression inputs (not plain objects) +// Verifies that Concat works with CwiseBinaryOp, CwiseUnaryOp, etc. +// ============================================================================ +template <typename MatrixType> +void test_concat_with_expressions(const MatrixType& m) { + typedef typename MatrixType::Scalar Scalar; + typedef Matrix<Scalar, Dynamic, Dynamic> MatrixX; + + Index rows = m.rows(); + Index cols = m.cols(); + + MatrixType a = MatrixType::Random(rows, cols); + MatrixType b = MatrixType::Random(rows, cols); + + // Concat with scalar-multiply expressions + { + MatrixX expected(2 * rows, cols); + expected.topRows(rows) = 2.0 * a; + expected.bottomRows(rows) = 3.0 * b; + VERIFY_IS_APPROX(expected, vcat(2.0 * a, 3.0 * b).eval()); + } + + // Concat with sum/difference expressions + { + MatrixX expected(rows, 2 * cols); + expected.leftCols(cols) = a + b; + expected.rightCols(cols) = a - b; + VERIFY_IS_APPROX(expected, hcat(a + b, a - b).eval()); + } + + // Concat with transpose expressions + { + MatrixX expected(cols, 2 * rows); + expected.leftCols(rows) = a.transpose(); + expected.rightCols(rows) = b.transpose(); + VERIFY_IS_APPROX(expected, hcat(a.transpose(), b.transpose()).eval()); + } + + // Concat with block expressions + { + if (rows >= 2 && cols >= 2) { + auto block_a = a.topLeftCorner(rows - 1, cols); + auto block_b = b.bottomRightCorner(1, cols); + MatrixX expected(rows, cols); + expected.topRows(rows - 1) = block_a; + expected.bottomRows(1) = block_b; + VERIFY_IS_APPROX(expected, vcat(block_a, block_b).eval()); + } + } +} + +// ============================================================================ +// Test 3: Fixed-size concat with compile-time dimension propagation +// ============================================================================ +template <typename Scalar, int LhsRows, int CommonCols, int RhsRows> +void test_vcat_fixed() { + typedef Matrix<Scalar, LhsRows, CommonCols> LhsType; + typedef Matrix<Scalar, RhsRows, CommonCols> RhsType; + + LhsType a = LhsType::Random(); + RhsType b = RhsType::Random(); + + auto vc = vcat(a, b); + VERIFY_IS_EQUAL(vc.rows(), a.rows() + b.rows()); + VERIFY_IS_EQUAL(vc.cols(), a.cols()); + + typedef Concat<Vertical, LhsType, RhsType> VConcatType; + VERIFY((int(VConcatType::RowsAtCompileTime) == LhsRows + RhsRows)); + VERIFY((int(VConcatType::ColsAtCompileTime) == CommonCols)); + VERIFY((int(VConcatType::MaxRowsAtCompileTime) == LhsRows + RhsRows)); + VERIFY((int(VConcatType::MaxColsAtCompileTime) == CommonCols)); + + typedef Matrix<Scalar, Dynamic, Dynamic> MatrixX; + MatrixX expected(a.rows() + b.rows(), a.cols()); + expected.topRows(a.rows()) = a; + expected.bottomRows(b.rows()) = b; + VERIFY_IS_APPROX(expected, MatrixX(vc)); +} + +template <typename Scalar, int CommonRows, int LhsCols, int RhsCols> +void test_hcat_fixed() { + typedef Matrix<Scalar, CommonRows, LhsCols> LhsType; + typedef Matrix<Scalar, CommonRows, RhsCols> RhsType; + + LhsType a = LhsType::Random(); + RhsType b = RhsType::Random(); + + auto hc = hcat(a, b); + VERIFY_IS_EQUAL(hc.rows(), a.rows()); + VERIFY_IS_EQUAL(hc.cols(), a.cols() + b.cols()); + + typedef Concat<Horizontal, LhsType, RhsType> HConcatType; + VERIFY((int(HConcatType::RowsAtCompileTime) == CommonRows)); + VERIFY((int(HConcatType::ColsAtCompileTime) == LhsCols + RhsCols)); + VERIFY((int(HConcatType::MaxRowsAtCompileTime) == CommonRows)); + VERIFY((int(HConcatType::MaxColsAtCompileTime) == LhsCols + RhsCols)); + + typedef Matrix<Scalar, Dynamic, Dynamic> MatrixX; + MatrixX expected(a.rows(), a.cols() + b.cols()); + expected.leftCols(a.cols()) = a; + expected.rightCols(b.cols()) = b; + VERIFY_IS_APPROX(expected, MatrixX(hc)); +} + +// ============================================================================ +// Test 4: Mixed fixed/dynamic dimension concat +// Verifies that compile-time dimensions become Dynamic when appropriate +// ============================================================================ +template <typename Scalar> +void test_concat_mixed_fixed_dynamic() { + typedef Matrix<Scalar, 3, 3> Fixed33; + typedef Matrix<Scalar, Dynamic, Dynamic> MatrixX; + typedef Matrix<Scalar, Dynamic, 3> MatrixX3; + + Fixed33 a = Fixed33::Random(); + MatrixX b = MatrixX::Random(4, 3); + + // Vertical: fixed + dynamic => dynamic rows, but cols known at compile time + { + auto vc = vcat(a, b); + VERIFY_IS_EQUAL(vc.rows(), 7); + VERIFY_IS_EQUAL(vc.cols(), 3); + + // Compile-time: rows should be Dynamic (fixed + dynamic = dynamic) + typedef decltype(vc) VCType; + VERIFY((int(VCType::RowsAtCompileTime) == Dynamic)); + // Cols: fixed 3 and dynamic => 3 (size_prefer_fixed picks the fixed one) + VERIFY((int(VCType::ColsAtCompileTime) == 3)); + + MatrixX expected(7, 3); + expected.topRows(3) = a; + expected.bottomRows(4) = b; + VERIFY_IS_APPROX(expected, MatrixX(vc)); + } + + // Vertical: fixed + MatrixX3 (cols both = 3) => dynamic rows + { + MatrixX3 c = MatrixX3::Random(5, 3); + auto vc = vcat(a, c); + VERIFY_IS_EQUAL(vc.rows(), 8); + VERIFY_IS_EQUAL(vc.cols(), 3); + + typedef decltype(vc) VCType; + VERIFY((int(VCType::RowsAtCompileTime) == Dynamic)); + + MatrixX expected(8, 3); + expected.topRows(3) = a; + expected.bottomRows(5) = c; + VERIFY_IS_APPROX(expected, MatrixX(vc)); + } + + // Horizontal: fixed + dynamic + { + MatrixX d = MatrixX::Random(3, 5); + auto hc = hcat(a, d); + VERIFY_IS_EQUAL(hc.rows(), 3); + VERIFY_IS_EQUAL(hc.cols(), 8); + + typedef decltype(hc) HCType; + VERIFY((int(HCType::ColsAtCompileTime) == Dynamic)); + + MatrixX expected(3, 8); + expected.leftCols(3) = a; + expected.rightCols(5) = d; + VERIFY_IS_APPROX(expected, MatrixX(hc)); + } + + // Dynamic + fixed (reversed order) + { + MatrixX3 e = MatrixX3::Random(5, 3); + auto vc = vcat(e, a); + VERIFY_IS_EQUAL(vc.rows(), 8); + VERIFY_IS_EQUAL(vc.cols(), 3); + + MatrixX expected(8, 3); + expected.topRows(5) = e; + expected.bottomRows(3) = a; + VERIFY_IS_APPROX(expected, MatrixX(vc)); + } +} + +// ============================================================================ +// Test 5: Rvalue temporary lifetime safety +// Verifies ref_selector properly copies temporaries instead of dangling refs +// ============================================================================ +template <typename Scalar> +void test_concat_rvalue_temporaries() { + typedef Matrix<Scalar, Dynamic, Dynamic> MatrixX; + + MatrixX a = MatrixX::Random(3, 3); + MatrixX b = MatrixX::Random(3, 3); + + // Test 5a: Product temporaries (Product expression creates temporaries) + { + MatrixX result = vcat(a * MatrixX::Identity(3, 3), b * MatrixX::Identity(3, 3)); + MatrixX expected(6, 3); + expected.topRows(3) = a; + expected.bottomRows(3) = b; + VERIFY_IS_APPROX(expected, result); + } + + // Test 5b: Additive temporaries + { + MatrixX result = hcat(a + MatrixX::Zero(3, 3), b + MatrixX::Zero(3, 3)); + MatrixX expected(3, 6); + expected.leftCols(3) = a; + expected.rightCols(3) = b; + VERIFY_IS_APPROX(expected, result); + } + + // Test 5c: Store expression with auto and evaluate later + // This is the critical rvalue lifetime test: the expression must + // keep the temporary alive until evaluation + { + auto expr = vcat(a * MatrixX::Identity(3, 3), b * MatrixX::Identity(3, 3)); + // At this point, the temporaries from a*I and b*I must still be valid + MatrixX result = expr; + MatrixX expected(6, 3); + expected.topRows(3) = a; + expected.bottomRows(3) = b; + VERIFY_IS_APPROX(expected, result); + } + + // Test 5d: Auto with hcat + { + auto expr = hcat(a + b, a - b); + MatrixX result = expr; + MatrixX expected(3, 6); + expected.leftCols(3) = a + b; + expected.rightCols(3) = a - b; + VERIFY_IS_APPROX(expected, result); + } + + // Test 5e: Nested auto expressions (chained concat of temporaries) + { + MatrixX c = MatrixX::Random(3, 3); + auto inner = vcat(a + b, c * MatrixX::Identity(3, 3)); + auto outer = vcat(inner, a - b); + MatrixX result = outer; + VERIFY_IS_EQUAL(result.rows(), 9); + VERIFY_IS_EQUAL(result.cols(), 3); + MatrixX expected(9, 3); + expected.topRows(3) = a + b; + expected.middleRows(3, 3) = c; + expected.bottomRows(3) = a - b; + VERIFY_IS_APPROX(expected, result); + } +} + +// ============================================================================ +// Test 6: Chained concatenation (binary API, no variadic needed) +// ============================================================================ +template <typename Scalar> +void test_concat_chained() { + typedef Matrix<Scalar, Dynamic, Dynamic> MatrixX; + + MatrixX a = MatrixX::Random(2, 3); + MatrixX b = MatrixX::Random(2, 3); + MatrixX c = MatrixX::Random(2, 3); + + // Chain: vcat(vcat(a, b), c) should be 6x3 + { + MatrixX result = vcat(vcat(a, b), c); + VERIFY_IS_EQUAL(result.rows(), 6); + VERIFY_IS_EQUAL(result.cols(), 3); + + MatrixX expected(6, 3); + expected.topRows(2) = a; + expected.middleRows(2, 2) = b; + expected.bottomRows(2) = c; + VERIFY_IS_APPROX(expected, result); + } + + // Chain: hcat(hcat(a, b), c) should be 2x9 + { + MatrixX result = hcat(hcat(a, b), c); + VERIFY_IS_EQUAL(result.rows(), 2); + VERIFY_IS_EQUAL(result.cols(), 9); + + MatrixX expected(2, 9); + expected.leftCols(3) = a; + expected.middleCols(3, 3) = b; + expected.rightCols(3) = c; + VERIFY_IS_APPROX(expected, result); + } + + // Mixed chain: vcat then hcat + { + MatrixX e = MatrixX::Random(4, 5); + MatrixX vert = vcat(a, b); // 4x3 + MatrixX result = hcat(vert, e); + VERIFY_IS_EQUAL(result.rows(), 4); + VERIFY_IS_EQUAL(result.cols(), 8); + + MatrixX expected(4, 8); + expected.leftCols(3) = vert; + expected.rightCols(5) = e; + VERIFY_IS_APPROX(expected, result); + } + + // 4-way chain + { + MatrixX d = MatrixX::Random(2, 3); + MatrixX result = vcat(vcat(vcat(a, b), c), d); + VERIFY_IS_EQUAL(result.rows(), 8); + VERIFY_IS_EQUAL(result.cols(), 3); + + MatrixX expected(8, 3); + expected.block(0, 0, 2, 3) = a; + expected.block(2, 0, 2, 3) = b; + expected.block(4, 0, 2, 3) = c; + expected.block(6, 0, 2, 3) = d; + VERIFY_IS_APPROX(expected, result); + } +} + +// ============================================================================ +// Test 7: Vector concatenation +// ============================================================================ +template <typename Scalar> +void test_concat_vectors() { + typedef Matrix<Scalar, Dynamic, 1> VectorX; + typedef Matrix<Scalar, 1, Dynamic> RowVectorX; + typedef Matrix<Scalar, 3, 1> Vector3; + typedef Matrix<Scalar, 4, 1> Vector4; + typedef Matrix<Scalar, 1, 3> RowVector3; + typedef Matrix<Scalar, 1, 4> RowVector4; + + // Column vectors: vcat + { + VectorX a = VectorX::Random(3); + VectorX b = VectorX::Random(5); + VectorX result = vcat(a, b); + VERIFY_IS_EQUAL(result.rows(), 8); + VERIFY_IS_EQUAL(result.cols(), 1); + VectorX expected(8); + expected.head(3) = a; + expected.tail(5) = b; + VERIFY_IS_APPROX(expected, result); + } + + // Row vectors: hcat + { + RowVectorX a = RowVectorX::Random(3); + RowVectorX b = RowVectorX::Random(5); + RowVectorX result = hcat(a, b); + VERIFY_IS_EQUAL(result.rows(), 1); + VERIFY_IS_EQUAL(result.cols(), 8); + RowVectorX expected(8); + expected.head(3) = a; + expected.tail(5) = b; + VERIFY_IS_APPROX(expected, result); + } + + // Fixed-size column vectors: compile-time dimensions + { + Vector3 a = Vector3::Random(); + Vector4 b = Vector4::Random(); + auto result = vcat(a, b); + + typedef Concat<Vertical, Vector3, Vector4> VConcatType; + VERIFY((int(VConcatType::RowsAtCompileTime) == 7)); + VERIFY((int(VConcatType::ColsAtCompileTime) == 1)); + + typedef Matrix<Scalar, Dynamic, 1> VecX; + VecX expected(7); + expected.head(3) = a; + expected.tail(4) = b; + VecX actual(result); + VERIFY_IS_APPROX(expected, actual); + } + + // Fixed-size row vectors: compile-time dimensions + { + RowVector3 a = RowVector3::Random(); + RowVector4 b = RowVector4::Random(); + auto result = hcat(a, b); + + typedef Concat<Horizontal, RowVector3, RowVector4> HConcatType; + VERIFY((int(HConcatType::RowsAtCompileTime) == 1)); + VERIFY((int(HConcatType::ColsAtCompileTime) == 7)); + + typedef Matrix<Scalar, 1, Dynamic> RowVecX; + RowVecX expected(7); + expected.head(3) = a; + expected.tail(4) = b; + RowVecX actual(result); + VERIFY_IS_APPROX(expected, actual); + } +} + +// ============================================================================ +// Test 8: Array concatenation (verifies XprKind propagation) +// ============================================================================ +template <typename Scalar> +void test_concat_array() { + typedef Array<Scalar, Dynamic, Dynamic> ArrayX; + typedef Array<Scalar, 2, 3> FixedArray23; + + // Dynamic arrays + { + ArrayX a = ArrayX::Random(3, 4); + ArrayX b = ArrayX::Random(3, 4); + + // vcat with arrays + ArrayX expected(6, 4); + expected.topRows(3) = a; + expected.bottomRows(3) = b; + ArrayX result = vcat(a, b); + VERIFY_IS_APPROX(expected, result); + } + + { + ArrayX a = ArrayX::Random(3, 4); + ArrayX b = ArrayX::Random(3, 4); + + // hcat with arrays + ArrayX expected(3, 8); + expected.leftCols(4) = a; + expected.rightCols(4) = b; + ArrayX result = hcat(a, b); + VERIFY_IS_APPROX(expected, result); + } + + // Fixed-size arrays + { + FixedArray23 a = FixedArray23::Random(); + FixedArray23 b = FixedArray23::Random(); + + typedef Concat<Vertical, FixedArray23, FixedArray23> VConcatType; + VERIFY((int(VConcatType::RowsAtCompileTime) == 4)); + VERIFY((int(VConcatType::ColsAtCompileTime) == 3)); + + ArrayX expected(4, 3); + expected.topRows(2) = a; + expected.bottomRows(2) = b; + ArrayX result = vcat(a, b); + VERIFY_IS_APPROX(expected, result); + } + + // Array expressions (coefficient-wise operations) + { + ArrayX a = ArrayX::Random(3, 4); + ArrayX b = ArrayX::Random(3, 4); + ArrayX result = vcat(a * b, a + b); + VERIFY_IS_EQUAL(result.rows(), 6); + VERIFY_IS_EQUAL(result.cols(), 4); + ArrayX expected(6, 4); + expected.topRows(3) = a * b; + expected.bottomRows(3) = a + b; + VERIFY_IS_APPROX(expected, result); + } +} + +// ============================================================================ +// Test 9: Coefficient-level access verification +// ============================================================================ +template <typename Scalar> +void test_concat_coeff_access() { + typedef Matrix<Scalar, Dynamic, Dynamic> MatrixX; + + MatrixX a = MatrixX::Random(3, 4); + MatrixX b = MatrixX::Random(2, 4); + + // Verify each coefficient of vcat + { + auto vc = vcat(a, b); + for (Index j = 0; j < 4; ++j) { + for (Index i = 0; i < 3; ++i) { + VERIFY_IS_APPROX(vc(i, j), a(i, j)); + } + for (Index i = 0; i < 2; ++i) { + VERIFY_IS_APPROX(vc(3 + i, j), b(i, j)); + } + } + } + + // Verify each coefficient of hcat + { + MatrixX c = MatrixX::Random(3, 2); + auto hc = hcat(a, c); + for (Index i = 0; i < 3; ++i) { + for (Index j = 0; j < 4; ++j) { + VERIFY_IS_APPROX(hc(i, j), a(i, j)); + } + for (Index j = 0; j < 2; ++j) { + VERIFY_IS_APPROX(hc(i, 4 + j), c(i, j)); + } + } + } +} + +// ============================================================================ +// Test 10: Expression used in further Eigen operations +// Verifies Concat integrates with Eigen's expression machinery +// ============================================================================ +template <typename Scalar> +void test_concat_in_expressions() { + typedef Matrix<Scalar, Dynamic, Dynamic> MatrixX; + typedef Matrix<Scalar, Dynamic, 1> VectorX; + + MatrixX a = MatrixX::Random(3, 3); + MatrixX b = MatrixX::Random(3, 3); + + // Use concat result in matrix multiplication + { + MatrixX wide = hcat(a, b); // 3x6 + VectorX v = VectorX::Random(6); + VectorX result = wide * v; + MatrixX wide_eval = wide; + VectorX expected = wide_eval * v; + VERIFY_IS_APPROX(expected, result); + } + + // Use concat result in sum + { + auto vc = vcat(a, b); + Scalar s = vc.sum(); + Scalar expected_sum = a.sum() + b.sum(); + VERIFY_IS_APPROX(expected_sum, s); + } + + // Use concat in norm computation + { + VectorX va = VectorX::Random(4); + VectorX vb = VectorX::Random(3); + auto vc = vcat(va, vb); + Scalar n2 = vc.squaredNorm(); + Scalar expected_n2 = va.squaredNorm() + vb.squaredNorm(); + VERIFY_IS_APPROX(expected_n2, n2); + } + + // Concat result assigned to a block + { + MatrixX target = MatrixX::Zero(6, 3); + target = vcat(a, b); + MatrixX expected(6, 3); + expected.topRows(3) = a; + expected.bottomRows(3) = b; + VERIFY_IS_APPROX(expected, target); + } +} + +// ============================================================================ +// Test 11: Single row / single column edge cases +// ============================================================================ +template <typename Scalar> +void test_concat_single_row_col() { + typedef Matrix<Scalar, Dynamic, Dynamic> MatrixX; + typedef Matrix<Scalar, 1, Dynamic> RowVectorX; + typedef Matrix<Scalar, Dynamic, 1> VectorX; + + // Single-row matrices: vcat + { + RowVectorX a = RowVectorX::Random(5); + RowVectorX b = RowVectorX::Random(5); + MatrixX result = vcat(a, b); + VERIFY_IS_EQUAL(result.rows(), 2); + VERIFY_IS_EQUAL(result.cols(), 5); + VERIFY_IS_APPROX(result.row(0), a); + VERIFY_IS_APPROX(result.row(1), b); + } + + // Single-column matrices: hcat + { + VectorX a = VectorX::Random(5); + VectorX b = VectorX::Random(5); + MatrixX result = hcat(a, b); + VERIFY_IS_EQUAL(result.rows(), 5); + VERIFY_IS_EQUAL(result.cols(), 2); + VERIFY_IS_APPROX(result.col(0), a); + VERIFY_IS_APPROX(result.col(1), b); + } + + // 1x1 matrices + { + MatrixX a = MatrixX::Random(1, 1); + MatrixX b = MatrixX::Random(1, 1); + MatrixX vr = vcat(a, b); + VERIFY_IS_EQUAL(vr.rows(), 2); + VERIFY_IS_EQUAL(vr.cols(), 1); + VERIFY_IS_APPROX(vr(0, 0), a(0, 0)); + VERIFY_IS_APPROX(vr(1, 0), b(0, 0)); + + MatrixX hr = hcat(a, b); + VERIFY_IS_EQUAL(hr.rows(), 1); + VERIFY_IS_EQUAL(hr.cols(), 2); + VERIFY_IS_APPROX(hr(0, 0), a(0, 0)); + VERIFY_IS_APPROX(hr(0, 1), b(0, 0)); + } +} + +// ============================================================================ +// Test 12: RowMajor matrices +// Verifies IsRowMajor flag logic in Concat traits +// ============================================================================ +template <typename Scalar> +void test_concat_row_major() { + typedef Matrix<Scalar, Dynamic, Dynamic, RowMajor> RowMajorMatrix; + typedef Matrix<Scalar, Dynamic, Dynamic> ColMajorMatrix; + typedef Matrix<Scalar, 3, 3, RowMajor> FixedRowMajor33; + + // RowMajor + RowMajor + { + RowMajorMatrix a = RowMajorMatrix::Random(3, 4); + RowMajorMatrix b = RowMajorMatrix::Random(2, 4); + ColMajorMatrix result = vcat(a, b); + VERIFY_IS_EQUAL(result.rows(), 5); + VERIFY_IS_EQUAL(result.cols(), 4); + ColMajorMatrix expected(5, 4); + expected.topRows(3) = a; + expected.bottomRows(2) = b; + VERIFY_IS_APPROX(expected, result); + } + + { + RowMajorMatrix a = RowMajorMatrix::Random(3, 4); + RowMajorMatrix b = RowMajorMatrix::Random(3, 2); + ColMajorMatrix result = hcat(a, b); + VERIFY_IS_EQUAL(result.rows(), 3); + VERIFY_IS_EQUAL(result.cols(), 6); + ColMajorMatrix expected(3, 6); + expected.leftCols(4) = a; + expected.rightCols(2) = b; + VERIFY_IS_APPROX(expected, result); + } + + // Mixed: RowMajor + ColMajor + { + RowMajorMatrix a = RowMajorMatrix::Random(3, 4); + ColMajorMatrix b = ColMajorMatrix::Random(2, 4); + ColMajorMatrix result = vcat(a, b); + VERIFY_IS_EQUAL(result.rows(), 5); + VERIFY_IS_EQUAL(result.cols(), 4); + ColMajorMatrix expected(5, 4); + expected.topRows(3) = a; + expected.bottomRows(2) = b; + VERIFY_IS_APPROX(expected, result); + } + + // Fixed-size RowMajor + { + FixedRowMajor33 a = FixedRowMajor33::Random(); + FixedRowMajor33 b = FixedRowMajor33::Random(); + ColMajorMatrix result = vcat(a, b); + VERIFY_IS_EQUAL(result.rows(), 6); + VERIFY_IS_EQUAL(result.cols(), 3); + ColMajorMatrix expected(6, 3); + expected.topRows(3) = a; + expected.bottomRows(3) = b; + VERIFY_IS_APPROX(expected, result); + } +} + +// ============================================================================ +// Test 13: Self-concatenation and aliasing +// ============================================================================ +template <typename Scalar> +void test_concat_self() { + typedef Matrix<Scalar, Dynamic, Dynamic> MatrixX; + typedef Matrix<Scalar, Dynamic, 1> VectorX; + + // Self-concat: vcat(a, a) + { + MatrixX a = MatrixX::Random(3, 4); + MatrixX result = vcat(a, a); + VERIFY_IS_EQUAL(result.rows(), 6); + VERIFY_IS_EQUAL(result.cols(), 4); + VERIFY_IS_APPROX(result.topRows(3), a); + VERIFY_IS_APPROX(result.bottomRows(3), a); + } + + // Self-concat: hcat(a, a) + { + MatrixX a = MatrixX::Random(3, 4); + MatrixX result = hcat(a, a); + VERIFY_IS_EQUAL(result.rows(), 3); + VERIFY_IS_EQUAL(result.cols(), 8); + VERIFY_IS_APPROX(result.leftCols(4), a); + VERIFY_IS_APPROX(result.rightCols(4), a); + } + + // Self-concat vector + { + VectorX v = VectorX::Random(5); + VectorX result = vcat(v, v); + VERIFY_IS_EQUAL(result.rows(), 10); + VERIFY_IS_APPROX(result.head(5), v); + VERIFY_IS_APPROX(result.tail(5), v); + } + + // Concat with self-expression: vcat(a, a.transpose().transpose()) + { + MatrixX a = MatrixX::Random(3, 4); + MatrixX result = vcat(a, a + MatrixX::Zero(3, 4)); + VERIFY_IS_EQUAL(result.rows(), 6); + VERIFY_IS_APPROX(result.topRows(3), a); + VERIFY_IS_APPROX(result.bottomRows(3), a); + } +} + +// ============================================================================ +// Test 14: Runtime assertion for mismatched dimensions +// ============================================================================ +template <typename Scalar> +void test_concat_mismatched_dimensions() { + typedef Matrix<Scalar, Dynamic, Dynamic> MatrixX; + + MatrixX a(3, 4); + MatrixX b(3, 5); + VERIFY_RAISES_ASSERT(vcat(a, b)); // cols don't match + + MatrixX c(3, 4); + MatrixX d(2, 4); + VERIFY_RAISES_ASSERT(hcat(c, d)); // rows don't match +} + +EIGEN_DECLARE_TEST(concat) { + for (int i = 0; i < g_repeat; i++) { + // Dynamic-size matrix concat + CALL_SUBTEST_1(test_concat_dynamic(MatrixXf(4, 3))); + CALL_SUBTEST_1(test_concat_dynamic(MatrixXd(5, 7))); + CALL_SUBTEST_1(test_concat_dynamic(MatrixXcf(3, 4))); + + // Expression inputs + CALL_SUBTEST_2(test_concat_with_expressions(MatrixXf(4, 3))); + CALL_SUBTEST_2(test_concat_with_expressions(MatrixXd(5, 7))); + + // Fixed-size concat with compile-time dimension checks + CALL_SUBTEST_3((test_vcat_fixed<float, 2, 3, 4>())); // vcat: 2x3 + 4x3 + CALL_SUBTEST_3((test_hcat_fixed<double, 3, 2, 5>())); // hcat: 3x2 + 3x5 + CALL_SUBTEST_3((test_vcat_fixed<float, 4, 4, 4>())); // square vcat: 4x4 + 4x4 + CALL_SUBTEST_3((test_hcat_fixed<float, 4, 4, 4>())); // square hcat: 4x4 + 4x4 + + // Mixed fixed/dynamic + CALL_SUBTEST_4(test_concat_mixed_fixed_dynamic<float>()); + CALL_SUBTEST_4(test_concat_mixed_fixed_dynamic<double>()); + + // Rvalue temporaries (lifetime safety) + CALL_SUBTEST_5(test_concat_rvalue_temporaries<float>()); + CALL_SUBTEST_5(test_concat_rvalue_temporaries<double>()); + + // Chained concat (binary API) + CALL_SUBTEST_6(test_concat_chained<float>()); + CALL_SUBTEST_6(test_concat_chained<double>()); + + // Vector concat + CALL_SUBTEST_7(test_concat_vectors<float>()); + CALL_SUBTEST_7(test_concat_vectors<double>()); + + // Array concat (XprKind propagation) + CALL_SUBTEST_8(test_concat_array<float>()); + CALL_SUBTEST_8(test_concat_array<double>()); + + // Coefficient-level access + CALL_SUBTEST_9(test_concat_coeff_access<float>()); + CALL_SUBTEST_9(test_concat_coeff_access<double>()); + + // Integration with other Eigen operations + CALL_SUBTEST_10(test_concat_in_expressions<float>()); + CALL_SUBTEST_10(test_concat_in_expressions<double>()); + + // Single row/column edge cases + CALL_SUBTEST_11(test_concat_single_row_col<float>()); + CALL_SUBTEST_11(test_concat_single_row_col<double>()); + + // RowMajor matrices + CALL_SUBTEST_12(test_concat_row_major<float>()); + CALL_SUBTEST_12(test_concat_row_major<double>()); + + // Self-concatenation and aliasing + CALL_SUBTEST_13(test_concat_self<float>()); + CALL_SUBTEST_13(test_concat_self<double>()); + + // Runtime assertion for mismatched dimensions + CALL_SUBTEST_14(test_concat_mismatched_dimensions<float>()); + CALL_SUBTEST_14(test_concat_mismatched_dimensions<double>()); + } +}