From a5034ce5ad99268c986ca8ebaed44c8c2a463df9 Mon Sep 17 00:00:00 2001
From: Mikhail Katliar <mikhail.katliar@protonmail.com>
Date: Fri, 25 Oct 2024 13:08:07 +0200
Subject: [PATCH] WIP: fixing StaticVectorPointerTest on neon64

---
 include/blast/math/dense/StaticVector.hpp     | 176 ++++++++++++
 include/blast/math/simd/arch/Neon64.hpp       |  17 +-
 include/blast/math/views/Row.hpp              | 261 ++++++++++++++++++
 .../math/dense/StaticVectorPointerTest.cpp    |  13 +-
 4 files changed, 449 insertions(+), 18 deletions(-)
 create mode 100644 include/blast/math/dense/StaticVector.hpp

diff --git a/include/blast/math/dense/StaticVector.hpp b/include/blast/math/dense/StaticVector.hpp
new file mode 100644
index 00000000..64c19b90
--- /dev/null
+++ b/include/blast/math/dense/StaticVector.hpp
@@ -0,0 +1,176 @@
+// Copyright (c) 2019-2024 Mikhail Katliar All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+#pragma once
+
+#include <blast/math/Forward.hpp>
+#include <blast/math/TransposeFlag.hpp>
+#include <blast/math/Simd.hpp>
+#include <blast/math/TypeTraits.hpp>
+#include <blast/system/CacheLine.hpp>
+#include <blast/util/NextMultiple.hpp>
+#include <blast/util/Types.hpp>
+
+#include <initializer_list>
+#include <type_traits>
+
+
+namespace blast
+{
+    /// @brief Vector with statically defined size.
+    ///
+    /// @tparam T element type of the vector
+    /// @tparam N number of elements
+    /// @tparam TF transpose flag
+    template <typename T, size_t N, TransposeFlag TF>
+    class StaticVector
+    {
+    public:
+        using ElementType = T;
+        static TransposeFlag constexpr transposeFlag = TF;
+
+
+        StaticVector() noexcept
+        {
+            // Initialize padding elements to 0 to prevent denorms in calculations.
+            // Denorms can significantly impair performance, see https://github.com/giaf/blasfeo/issues/103
+            std::fill_n(v_, capacity_, T {});
+        }
+
+
+        StaticVector(T const& v) noexcept
+        {
+            std::fill_n(v_, capacity_, v);
+        }
+
+
+        /**
+         * @brief Construct from an initializer list.
+         *
+         * \code
+         * StaticVector<double, 3, columnVector> v {1., 2., 3.};
+         * \endcode
+         *
+         * @param list list of vector elements. If @a list is shorter than @a N,
+         * the remaining vector elements will be 0. If @a list is longer than @a N,
+         * the extra elements of @a list will be ignored.
+         */
+        constexpr StaticVector(std::initializer_list<T> list)
+        {
+            fill(copy_n(list.begin(), std::min(list.size(), N), std::begin(v_)), std::end(v_), T {});
+        }
+
+
+        StaticVector& operator=(T val) noexcept
+        {
+            std::fill_n(v_, capacity_, val);
+
+            return *this;
+        }
+
+
+        constexpr T const& operator[](size_t i) const noexcept
+        {
+            assert(i < N);
+            return v_[i];
+        }
+
+
+        constexpr T& operator[](size_t i)
+        {
+            assert(i < N);
+            return v_[i];
+        }
+
+
+        static size_t constexpr size() noexcept
+        {
+            return N;
+        }
+
+
+        T * data() noexcept
+        {
+            return v_;
+        }
+
+
+        T const * data() const noexcept
+        {
+            return v_;
+        }
+
+
+        /**
+         * @brief Set all vector elements to 0
+         */
+        void reset() noexcept
+        {
+            std::fill_n(v_, capacity_, T {});
+        }
+
+
+    private:
+        static size_t constexpr capacity_ = nextMultiple(N, SimdSize_v<T>);
+
+        // Alignment of the data elements.
+        static size_t constexpr alignment_ = CACHE_LINE_SIZE;
+
+        // Aligned element storage.
+        alignas(alignment_) T v_[capacity_];
+
+    };
+
+
+    template <typename T, size_t N, TransposeFlag TF>
+    inline size_t constexpr size(StaticVector<T, N, TF> const& m) noexcept
+    {
+        return N;
+    }
+
+
+    template <typename T, size_t N, TransposeFlag TF>
+    inline constexpr T * data(StaticVector<T, N, TF>& m) noexcept
+    {
+        return m.data();
+    }
+
+
+    template <typename T, size_t N, TransposeFlag TF>
+    inline constexpr T const * data(StaticVector<T, N, TF> const& m) noexcept
+    {
+        return m.data();
+    }
+
+
+    template <typename T, size_t N, TransposeFlag TF>
+    inline void reset(StaticVector<T, N, TF>& m) noexcept
+    {
+        m.reset();
+    }
+
+
+    template <typename T, size_t N, TransposeFlag TF>
+    struct IsDenseVector<StaticVector<T, N, TF>> : std::true_type {};
+
+
+    template <typename T, size_t N, TransposeFlag TF>
+    struct IsStatic<StaticVector<T, N, TF>> : std::true_type {};
+
+
+    template <typename T, size_t N, TransposeFlag TF>
+    struct IsAligned<StaticVector<T, N, TF>> : std::integral_constant<bool, true> {};
+
+
+    template <typename T, size_t N, TransposeFlag TF>
+    struct IsPadded<StaticVector<T, N, TF>> : std::integral_constant<bool, true> {};
+
+
+    template <typename T, size_t N, TransposeFlag TF>
+    struct IsStaticallySpaced<StaticVector<T, N, TF>> : std::integral_constant<bool, true> {};
+
+
+    template <typename T, size_t N, TransposeFlag TF>
+    struct Spacing<StaticVector<T, N, TF>> : std::integral_constant<size_t, 1> {};
+}
diff --git a/include/blast/math/simd/arch/Neon64.hpp b/include/blast/math/simd/arch/Neon64.hpp
index 83c23cbc..6c2238f9 100644
--- a/include/blast/math/simd/arch/Neon64.hpp
+++ b/include/blast/math/simd/arch/Neon64.hpp
@@ -1,16 +1,7 @@
-// Copyright 2024 Mikhail Katliar
-//
-// Licensed under the Apache License, Version 2.0 (the "License");
-// you may not use this file except in compliance with the License.
-// You may obtain a copy of the License at
-//
-//     https://www.apache.org/licenses/LICENSE-2.0
-//
-// Unless required by applicable law or agreed to in writing, software
-// distributed under the License is distributed on an "AS IS" BASIS,
-// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
-// See the License for the specific language governing permissions and
-// limitations under the License.
+// Copyright (c) 2019-2024 Mikhail Katliar All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
 #pragma once
 
 #include <xsimd/xsimd.hpp>
diff --git a/include/blast/math/views/Row.hpp b/include/blast/math/views/Row.hpp
index b9afee98..9cd61fc6 100644
--- a/include/blast/math/views/Row.hpp
+++ b/include/blast/math/views/Row.hpp
@@ -6,3 +6,264 @@
 
 #include <blast/math/views/row/BaseTemplate.hpp>
 #include <blast/math/views/row/Panel.hpp>
+#include <blast/math/TypeTraits.hpp>
+
+
+namespace blast
+{
+    /**
+     * @brief Row view of a matrix
+     *
+     * NOTE: this implementation is not optimized!
+     *
+     * @tparam MT viewed matrix type
+     */
+    template <Matrix MT>
+    class Row
+    {
+    public:
+        static TransposeFlag constexpr transposeFlag = rowVector;
+
+        using ViewedType = MT;
+        using ElementType = ElementType_t<MT>;
+
+        //! Reference to a constant submatrix value.
+        using ConstReference = ElementType const&;
+
+        //! Reference to a non-constant submatrix value.
+        using Reference = std::conditional_t<IsConst_v<MT>, ConstReference, ElementType&>;
+
+        //! Pointer to a constant submatrix value.
+        using ConstPointer = ElementType const *;
+
+        //! Pointer to a non-constant submatrix value.
+        using Pointer = std::conditional_t<IsConst_v<MT>, ConstPointer, ElementType *>;
+
+
+        /**
+         * @brief Constructor
+         *
+         * @param matrix the matrix
+         * @param i row index
+         * @param j start of the row
+         * @param n row length
+         */
+        explicit inline constexpr Row(MT& matrix, size_t i, size_t j, size_t n)
+        :   matrix_ {matrix}
+        ,   i_ {i}
+        ,   j_ {j}
+        ,   n_ {n}
+        ,   data_ {&matrix(i, j)}
+        {
+            BLAST_USER_ASSERT(i_ < rows(matrix_), "Invalid row index");
+            BLAST_USER_ASSERT(j_ < columns(matrix_), "Invalid column index");
+            BLAST_USER_ASSERT(j_ + n_ <= columns(matrix_), "Invalid row length");
+        }
+
+
+        Row(Row const&) = default;
+
+
+        /**
+         * @brief Vector assignment
+         *
+         * Copies elements from the right-hand side expression
+         *
+         * @tparam VT type of right-hand side vector
+         *
+         * @param rhs
+         *
+         * @return reference to *this
+         */
+        template <Vector VT>
+        Row& operator=(VT const& rhs)
+        {
+            assign(*this, rhs);
+            return *this;
+        }
+
+
+        Reference operator[](size_t j) noexcept
+        {
+            BLAST_USER_ASSERT(j_ + j < columns(matrix_), "Invalid vector access index");
+
+            return matrix_(i_, j_ + j);
+        }
+
+
+        ConstReference operator[](size_t j) const
+        {
+            BLAST_USER_ASSERT(j_ + j < columns(matrix_), "Invalid vector access index");
+
+            return const_cast<MT const&>(matrix_)(i_, j_ + j);
+        }
+
+
+        friend Pointer data(Row& m) noexcept
+        {
+            return m.data_;
+        }
+
+
+        friend ConstPointer data(Row const& m) noexcept
+        {
+            return m.data_;
+        }
+
+
+        friend size_t constexpr size() const noexcept
+        {
+            return n_;
+        }
+
+
+        /**
+         * @brief Subvector of a row
+         *
+         * @param row the original row
+         * @param j start index of the subvector within @a row
+         * @param n length of the subvector
+         */
+        friend auto subvector(Row& row, size_t j, size_t n) noexcept
+        {
+            return Row<MT> {row.matrix_, row.i_, row.j_ + j, n};
+        }
+
+
+        /**
+         * @brief Subvector of a const row
+         *
+         * @param row the original row
+         * @param j start index of the subvector within @a row
+         * @param n length of the subvector
+         */
+        friend auto subvector(Row const& row, size_t j, size_t n) noexcept
+        {
+            return Row<MT const> {const_cast<MT const&>(row.matrix_), row.i_, row.j_ + j, n};
+        }
+
+
+    private:
+        ViewedType& matrix_;        //!< The matrix containing the submatrix.
+        size_t const i_;
+        size_t const j_;
+        size_t const n_;
+        Pointer const data_;
+    };
+
+
+    /**
+     * @brief Specialization for @a Row class
+     */
+    template <typename MT>
+    struct IsAligned<Row<MT>> : std::integral_constant<bool, IsAligned_v<MT> && IsPadded_v<MT>> {};
+
+
+    /**
+     * @brief Specialization for @a Row class
+     */
+    template <typename MT>
+    struct IsPadded<Row<MT>> : std::integral_constant<bool, IsPadded_v<MT> && StorageOrder_v<MT> == rowMajor> {};
+
+
+    /**
+     * @brief Specialization for @a Row class
+     */
+    template <typename MT>
+    struct IsStatic<Row<MT>> : IsStatic<MT> {};
+
+
+    /**
+     * @brief Specialization for @a Row class
+     */
+    template <typename MT>
+    struct IsDenseVector<Row<MT>> : IsDenseMatrix<MT> {};
+
+
+    /**
+     * @brief Specialization for @a Row class
+     */
+    template <typename MT>
+    struct IsView<Row<MT>> : std::integral_constant<bool, true> {};
+
+
+    /**
+     * @brief Specialization for rows of dense (non-panel) matrices
+     */
+    template <typename MT>
+    requires IsDenseMatrix_v<MT>
+    struct IsStaticallySpaced<Row<MT>> : std::integral_constant<bool, IsStatic_v<MT> || StorageOrder_v<MT> == rowMajor> {};
+
+
+    /**
+     * @brief Specialization for rows of panel matrices
+     */
+    template <typename MT>
+    requires IsPanelMatrix_v<MT>
+    struct IsStaticallySpaced<Row<MT>> : std::integral_constant<bool, StorageOrder_v<MT> == columnMajor> {};
+
+
+    /**
+     * @brief Specialization for rows of dense (non-panel) matrices
+     */
+    template <typename MT>
+    requires IsDenseMatrix_v<MT> && IsStatic_v<MT>
+    struct Spacing<Row<MT>> : std::integral_constant<bool, StorageOrder_v<MT> == rowMajor ? 1 : Spacing_v<MT>> {};
+
+
+    /**
+     * @brief Specialization for rows of panel matrices
+     */
+    template <typename MT>
+    requires IsPanelMatrix_v<MT> && (StorageOrder_v<MT> == columnMajor)
+    struct Spacing<Row<MT>> : std::integral_constant<bool, SimdSize_v<ElementType_t<MT>>> {};
+
+
+    /**
+     * @brief Full row of a matrix
+     *
+     * @tparam MT type of the matrix containing the row
+     *
+     * @param matrix matrix containing the row
+     * @param row row index
+     *
+     * @return submatrix of @a matrix
+     */
+    template <typename MT>
+    inline auto row(MT& matrix, size_t row)
+    {
+        return Row<MT> {matrix, row, 0, columns(matrix)};
+    }
+
+
+    /**
+     * @brief Partial row of a matrix
+     *
+     * @tparam MT type of the matrix containing the row
+     *
+     * @param matrix matrix containing the row
+     * @param row row index
+     * @param column row start
+     * @param n row length
+     *
+     * @return submatrix of @a matrix
+     */
+    template <typename MT>
+    inline auto row(MT& matrix, size_t row, size_t column, size_t n)
+    {
+        return Row<MT> {matrix, row, column, n};
+    }
+
+
+    /**
+     * @brief Set all elements of a row to their default value (0).
+     *
+     * @param matrix submatrix to set to 0.
+     */
+    template <typename MT>
+    inline void reset(Row<MT>& row) noexcept
+    {
+        for (size_t i = 0; i < size(row); ++i)
+            row[i] = 0;
+    }
+}
diff --git a/test/blast/math/dense/StaticVectorPointerTest.cpp b/test/blast/math/dense/StaticVectorPointerTest.cpp
index a8050c8f..b737aa9b 100644
--- a/test/blast/math/dense/StaticVectorPointerTest.cpp
+++ b/test/blast/math/dense/StaticVectorPointerTest.cpp
@@ -3,7 +3,10 @@
 // license that can be found in the LICENSE file.
 
 #include <blast/math/Vector.hpp>
-#include <blast/blaze/Math.hpp>
+#include <blast/math/Matrix.hpp>
+#include <blast/math/dense/StaticVector.hpp>
+#include <blast/math/dense/StaticMatrix.hpp>
+#include <blast/math/views/Row.hpp>
 
 #include <test/Testing.hpp>
 
@@ -18,7 +21,7 @@ namespace blast :: testing
         using Real = Scalar;
 
 
-        template <bool TF>
+        template <TransposeFlag TF>
         void testSpacingImpl()
         {
             StaticVector<Real, 3, TF> v;
@@ -27,7 +30,7 @@ namespace blast :: testing
         }
 
 
-        template <bool TF>
+        template <TransposeFlag TF>
         void testGetImpl()
         {
             StaticVector<Real, 3, TF> v;
@@ -37,7 +40,7 @@ namespace blast :: testing
         }
 
 
-        template <bool TF>
+        template <TransposeFlag TF>
         void testOffsetImpl()
         {
             StaticVector<Real, 5, TF> v;
@@ -56,7 +59,7 @@ namespace blast :: testing
             StaticMatrix<Real, 5, 5, SO> A;
 
             size_t constexpr i = 1, j = 2;
-            auto p = ptr<unaligned>(blaze::subvector<j, columns(A) - j>(blaze::row<i>(A)), 0);
+            auto p = ptr<unaligned>(subvector(row(A, i), j, columns(A) - j), 0);
             ASSERT_EQ(p.get(), &A(i, j));
             ASSERT_EQ(p.spacing(), SO == columnMajor ? A.spacing() : 1);
         }