eckit::linalg::SparseMatrix::Allocator StandardContainerAllocator

Author: pmaciel

src/eckit/linalg/CMakeLists.txt

@@ -24,6 +24,8 @@ list( APPEND eckit_la_srcs
       allocator/NonOwningAllocator.h
       allocator/StandardAllocator.cc
       allocator/StandardAllocator.h
+      allocator/StandardContainerAllocator.cc
+      allocator/StandardContainerAllocator.h
       dense/LinearAlgebraGeneric.cc
       dense/LinearAlgebraGeneric.h
       sparse/LinearAlgebraGeneric.cc

src/eckit/linalg/allocator/StandardContainerAllocator.cc

@@ -0,0 +1,78 @@
+/*
+ * (C) Copyright 2024- ECMWF.
+ *
+ * This software is licensed under the terms of the Apache Licence Version 2.0
+ * which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
+ * In applying this licence, ECMWF does not waive the privileges and immunities
+ * granted to it by virtue of its status as an intergovernmental organisation
+ * nor does it submit to any jurisdiction.
+ */
+
+
+#include "eckit/linalg/allocator/StandardContainerAllocator.h"
+
+#include <numeric>
+#include <ostream>
+
+#include "eckit/exception/Exceptions.h"
+
+
+namespace eckit::linalg::allocator {
+
+
+StandardContainerAllocator::StandardContainerAllocator(Size Nr, Size Nc, const container_type& cont) :
+    Nr_(Nr), Nc_(Nc), nnz_([Nc](const auto& cont) {
+        return std::accumulate(cont.begin(), cont.end(), 0, [Nc](const auto& acc, const auto& row) {
+            ASSERT(row.empty() || row.rbegin()->first < Nc);
+            return acc + row.size();
+        });
+    }(cont)) {
+    ASSERT(Nr == cont.size());
+
+    ia_.reserve(Nr_ + 1);
+    ja_.reserve(nnz_);
+    a_.reserve(nnz_);
+
+    ia_.emplace_back(0);
+    for (const auto& row : cont) {
+        ia_.emplace_back(ia_.back() + row.size());
+        for (const auto& [col, val] : row) {
+            ja_.emplace_back(static_cast<Index>(col));
+            a_.emplace_back(val);
+        }
+    }
+}
+
+
+SparseMatrix::Layout StandardContainerAllocator::StandardContainerAllocator::allocate(SparseMatrix::Shape& shape) {
+    shape.size_ = nnz_;
+    shape.rows_ = Nr_;
+    shape.cols_ = Nc_;
+
+    SparseMatrix::Layout layout;
+    layout.outer_ = ia_.data();
+    layout.inner_ = ja_.data();
+    layout.data_  = a_.data();
+
+    return layout;
+}
+
+
+void StandardContainerAllocator::deallocate(SparseMatrix::Layout, SparseMatrix::Shape) {
+    ia_.clear();
+    ja_.clear();
+    a_.clear();
+}
+
+
+bool StandardContainerAllocator::inSharedMemory() const {
+    return false;
+}
+
+
+void StandardContainerAllocator::print(std::ostream& out) const {
+    out << "StandardContainerAllocator[Nr=" << Nr_ << ",Nc=" << Nc_ << ",nnz=" << nnz_ << "]";
+}
+
+
+}  // namespace eckit::linalg::allocator

src/eckit/linalg/allocator/StandardContainerAllocator.h

@@ -0,0 +1,50 @@
+/*
+ * (C) Copyright 2024- ECMWF.
+ *
+ * This software is licensed under the terms of the Apache Licence Version 2.0
+ * which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
+ * In applying this licence, ECMWF does not waive the privileges and immunities
+ * granted to it by virtue of its status as an intergovernmental organisation
+ * nor does it submit to any jurisdiction.
+ */
+
+
+#pragma once
+
+#include <map>
+#include <vector>
+
+#include "eckit/linalg/SparseMatrix.h"
+
+
+namespace eckit::linalg::allocator {
+
+
+/**
+ * @brief In-place allocator for sparse matrices, directly mapping supporting arrays including from another
+ * SparseMatrix. It is able to provide a "view".
+ */
+class StandardContainerAllocator : public SparseMatrix::Allocator {
+public:
+    using container_type = std::vector<std::map<Size, Scalar>>;
+
+    StandardContainerAllocator(Size Nr, Size Nc, const container_type&);
+
+    SparseMatrix::Layout allocate(SparseMatrix::Shape&) override;
+
+    void deallocate(SparseMatrix::Layout, SparseMatrix::Shape) override;
+    bool inSharedMemory() const override;
+    void print(std::ostream&) const override;
+
+private:
+    const Size Nr_;
+    const Size Nc_;
+    const Size nnz_;
+
+    std::vector<SparseMatrix::UIndex> ia_;
+    std::vector<Index> ja_;
+    std::vector<Scalar> a_;
+};
+
+
+}  // namespace eckit::linalg::allocator

tests/linalg/test_la_sparse.cc

@@ -12,6 +12,7 @@
 #include "util.h"
 
 #include "eckit/linalg/allocator/NonOwningAllocator.h"
+#include "eckit/linalg/allocator/StandardContainerAllocator.h"
 
 using namespace eckit::linalg;
 
@@ -212,6 +213,16 @@ CASE("eckit la sparse") {
         EXPECT(G.A.inner() == GA.inner());
         EXPECT(G.A.data() == GA.data());
     }
+
+    SECTION("containers allocator") {
+        SparseMatrix FA(new allocator::StandardContainerAllocator(3, 3, {{{0, 2.}, {2, -3.}}, {{1, 2.}}, {{2, 2.}}}));
+
+        EXPECT(equal_sparse_matrix(FA, F.A.outer(), F.A.inner(), F.A.data()));
+
+        SparseMatrix GA(new allocator::StandardContainerAllocator(2, 3, {{{0, 1.}, {2, 2.}}, {{0, 3.}, {1, 4.}}}));
+
+        EXPECT(equal_sparse_matrix(GA, GA.outer(), G.A.inner(), G.A.data()));
+    }
 }
 
 //----------------------------------------------------------------------------------------------------------------------