Detailed Description

template<std::size_t Dim>
class ProcessLib::Assembly::MatrixElementCache< Dim >

Definition at line 65 of file MatrixElementCache.h.

#include <MatrixElementCache.h>

Collaboration diagram for ProcessLib::Assembly::MatrixElementCache< Dim >:

Public Member Functions
	MatrixElementCache (MatOrVec &mat_or_vec, Stats &stats, const int num_threads)
void	add (std::vector< double > const &local_data, std::vector< GlobalIndexType > const &indices)
	~MatrixElementCache ()

Private Types
using	MatOrVec = std::conditional_t<Dim == 1, GlobalVector, GlobalMatrix>

Private Member Functions
void	addToCache (std::vector< double > const &values, std::vector< GlobalIndexType > const &indices)
void	addToCacheImpl (std::vector< double > const &values, std::vector< GlobalIndexType > const &indices, std::integral_constant< std::size_t, 1 >)
void	addToCacheImpl (std::vector< double > const &values, std::vector< GlobalIndexType > const &indices, std::integral_constant< std::size_t, 2 >)
void	ensureEnoughSpace (std::size_t const space_needed)
OGS_NO_INLINE void	ensureEnoughSpaceSlow (std::size_t const space_needed)
void	addCacheToGlobal ()

Private Attributes
std::vector< MatrixElementCacheEntry< Dim > >	cache_
MatOrVec &	mat_or_vec_
Stats &	stats_
int	num_threads_

Static Private Attributes
static constexpr std::size_t	cache_capacity = 1'000'000

Member Typedef Documentation

◆ MatOrVec

template<std::size_t Dim>

using ProcessLib::Assembly::MatrixElementCache< Dim >::MatOrVec = std::conditional_t<Dim == 1, GlobalVector, GlobalMatrix>

private

Definition at line 68 of file MatrixElementCache.h.

Constructor & Destructor Documentation

◆ MatrixElementCache()

template<std::size_t Dim>

ProcessLib::Assembly::MatrixElementCache< Dim >::MatrixElementCache	(	MatOrVec &	mat_or_vec,
		Stats &	stats,
		const int	num_threads )

inline

Definition at line 73 of file MatrixElementCache.h.

        : mat_or_vec_(mat_or_vec), stats_(stats), num_threads_(num_threads)
    {
        cache_.reserve(cache_capacity);
    }

References cache_, cache_capacity, mat_or_vec_, num_threads_, and stats_.

◆ ~MatrixElementCache()

template<std::size_t Dim>

ProcessLib::Assembly::MatrixElementCache< Dim >::~MatrixElementCache ( )

inline

Definition at line 108 of file MatrixElementCache.h.

108{ addCacheToGlobal(); }

ProcessLib::Assembly::MatrixElementCache::addCacheToGlobal

void addCacheToGlobal()

Definition MatrixElementCache.h:227

References addCacheToGlobal().

Member Function Documentation

◆ add()

template<std::size_t Dim>

void ProcessLib::Assembly::MatrixElementCache< Dim >::add	(	std::vector< double > const &	local_data,
		std::vector< GlobalIndexType > const &	indices )

inline

Definition at line 80 of file MatrixElementCache.h.

    {
        if (local_data.empty())
        {
            return;
        }
 
        if (num_threads_ == 1)
        {
            if constexpr (Dim == 2)
            {
                auto const num_r_c = indices.size();
                mat_or_vec_.add(
                    MathLib::RowColumnIndices<GlobalIndexType>{indices,
                                                               indices},
                    MathLib::toMatrix(local_data, num_r_c, num_r_c));
            }
            else
            {
                mat_or_vec_.add(indices, local_data);
            }
            return;
        }
 
        addToCache(local_data, indices);
    }

References addToCache(), mat_or_vec_, num_threads_, and MathLib::toMatrix().

◆ addCacheToGlobal()

template<std::size_t Dim>

void ProcessLib::Assembly::MatrixElementCache< Dim >::addCacheToGlobal ( )

inlineprivate

Definition at line 227 of file MatrixElementCache.h.

    {
        if (cache_.empty())
        {
            return;
        }
 
#pragma omp critical
        {
            if constexpr (Dim == 2)
            {
                auto const n_cols = mat_or_vec_.getNumberOfColumns();
 
                // TODO would be more efficient if our global matrix and vector
                // implementations supported batch addition of matrix elements
                // with arbitrary indices (not restricted to (n x m) shaped
                // submatrices).
                for (auto const [rc, value] : cache_)
                {
                    auto const [r, c] = rc;
 
                    auto const c_no_ghost =
                        detail::transformToNonGhostIndex(c, n_cols);
 
                    mat_or_vec_.add(r, c_no_ghost, value);
                }
            }
            else
            {
                // TODO batch addition would be more efficient. That needs the
                // refactoring of the matrix element cache.
                for (auto const [r, value] : cache_)
                {
                    mat_or_vec_.add(r.front(), value);
                }
            }
        }
 
        stats_.count_global += cache_.size();
        cache_.clear();
    }

References cache_, mat_or_vec_, stats_, and ProcessLib::Assembly::detail::transformToNonGhostIndex().

Referenced by ~MatrixElementCache(), and ensureEnoughSpaceSlow().

◆ addToCache()

template<std::size_t Dim>

void ProcessLib::Assembly::MatrixElementCache< Dim >::addToCache	(	std::vector< double > const &	values,
		std::vector< GlobalIndexType > const &	indices )

inlineprivate

Definition at line 111 of file MatrixElementCache.h.

    {
        ensureEnoughSpace(values.size());
 
        addToCacheImpl(values, indices,
                       std::integral_constant<std::size_t, Dim>{});
    }

References addToCacheImpl(), and ensureEnoughSpace().

Referenced by add().

◆ addToCacheImpl() [1/2]

template<std::size_t Dim>

void ProcessLib::Assembly::MatrixElementCache< Dim >::addToCacheImpl	(	std::vector< double > const &	values,
		std::vector< GlobalIndexType > const &	indices,
		std::integral_constant< std::size_t, 1 >	)

inlineprivate

Definition at line 121 of file MatrixElementCache.h.

    {
        auto const num_r_c = indices.size();
        auto const* const __restrict val = values.data();
        auto const* const __restrict idx = indices.data();
 
        // Count non-zeros first
        auto const nonzero_count = ranges::count_if(
            val, val + num_r_c, [](double v) { return v != 0.0; });
 
        if (nonzero_count == 0)
        {
            stats_.count += num_r_c;
            return;
        }
 
        auto entries =
            ranges::views::iota(0ul, num_r_c) |
            ranges::views::filter([val](std::size_t i)
                                  { return val[i] != 0.0; }) |
            ranges::views::transform(
                [val, idx](std::size_t i)
                { return MatrixElementCacheEntry<1>{{idx[i]}, val[i]}; });
 
        auto const old_size = cache_.size();
        cache_.resize(old_size + nonzero_count);
        ranges::copy(entries, cache_.begin() + old_size);
 
        stats_.count += num_r_c;
        stats_.count_nonzero += nonzero_count;
    }

References cache_, and stats_.

Referenced by addToCache().

◆ addToCacheImpl() [2/2]

template<std::size_t Dim>

void ProcessLib::Assembly::MatrixElementCache< Dim >::addToCacheImpl	(	std::vector< double > const &	values,
		std::vector< GlobalIndexType > const &	indices,
		std::integral_constant< std::size_t, 2 >	)

inlineprivate

Definition at line 156 of file MatrixElementCache.h.

    {
        auto const num_r_c = indices.size();
        auto const total_entries = num_r_c * num_r_c;
 
        // Store pointers to avoid repeated dereferencing.
        auto const* __restrict__ const idx = indices.data();
        auto const* __restrict__ const val = values.data();
 
        auto const old_size = cache_.size();
        cache_.resize(old_size + total_entries);
        auto* __restrict__ dest = cache_.data() + old_size;
 
        for (std::size_t r = 0; r < num_r_c; ++r)
        {
            auto const r_global = idx[r];
            for (std::size_t c = 0; c < num_r_c; ++c)
            {
                dest->indices[0] = r_global;
                dest->indices[1] = idx[c];
                dest->value = val[r * num_r_c + c];  // Row-major indexing
                ++dest;
            }
        }
 
        stats_.count += total_entries;
        stats_.count_nonzero += total_entries;
    }

References cache_, and stats_.

◆ ensureEnoughSpace()

template<std::size_t Dim>

void ProcessLib::Assembly::MatrixElementCache< Dim >::ensureEnoughSpace ( std::size_t const space_needed )

inlineprivate

Definition at line 187 of file MatrixElementCache.h.

    {
        // Fast path - inline this check
        if (cache_.size() + space_needed <= cache_.capacity()) [[likely]]
        {
            return;
        }
 
        // Slow path - keep in separate function to avoid code bloat
        ensureEnoughSpaceSlow(space_needed);
    }

References cache_, and ensureEnoughSpaceSlow().

Referenced by addToCache().

◆ ensureEnoughSpaceSlow()

template<std::size_t Dim>

OGS_NO_INLINE void ProcessLib::Assembly::MatrixElementCache< Dim >::ensureEnoughSpaceSlow ( std::size_t const space_needed )

inlineprivate

Definition at line 199 of file MatrixElementCache.h.

    {
        if (cache_.capacity() < cache_capacity)
        {
            cache_.reserve(cache_capacity);
        }
 
        auto const size_initial = cache_.size();
        auto const cap_initial = cache_.capacity();
 
        if (size_initial + space_needed <= cap_initial)
        {
            return;
        }
 
        addCacheToGlobal();
 
        // ensure again that there is enough capacity (corner case, if initial
        // capacity is too small because of whatever arcane reason)
        auto const size_after = cache_.size();
        auto const cap_after = cache_.capacity();
 
        if (size_after + space_needed > cap_after)
        {
            cache_.reserve(size_after + 2 * space_needed);
        }
    }