[PATCH v1 5/6] libstdc++: Implement layout_stride from mdspan.

Luc Grosheintz luc.grosheintz@gmail.com
Tue May 20 10:14:49 GMT 2025



On 5/20/25 11:40 AM, Tomasz Kaminski wrote:
> On Tue, May 20, 2025 at 11:20 AM Luc Grosheintz <luc.grosheintz@gmail.com>
> wrote:
> 
>>
>>
>> On 5/20/25 10:48 AM, Tomasz Kaminski wrote:
>>> On Tue, May 20, 2025 at 10:45 AM Luc Grosheintz <
>> luc.grosheintz@gmail.com>
>>> wrote:
>>>
>>>>
>>>>
>>>> On 5/20/25 10:24 AM, Tomasz Kaminski wrote:
>>>>> On Sun, May 18, 2025 at 10:16 PM Luc Grosheintz <
>>>> luc.grosheintz@gmail.com>
>>>>> wrote:
>>>>>
>>>>>> Implements the remaining parts of layout_left and layout_right; and
>> all
>>>>>> of layout_stride.
>>>>>>
>>>>>> libstdc++-v3/ChangeLog:
>>>>>>
>>>>>>            * include/std/mdspan(layout_stride): New class.
>>>>>>
>>>>>> Signed-off-by: Luc Grosheintz <luc.grosheintz@gmail.com>
>>>>>> ---
>>>>>>     libstdc++-v3/include/std/mdspan | 219
>> +++++++++++++++++++++++++++++++-
>>>>>>     1 file changed, 216 insertions(+), 3 deletions(-)
>>>>>>
>>>>>> diff --git a/libstdc++-v3/include/std/mdspan
>>>>>> b/libstdc++-v3/include/std/mdspan
>>>>>> index b1984eb2a33..31a38c736c2 100644
>>>>>> --- a/libstdc++-v3/include/std/mdspan
>>>>>> +++ b/libstdc++-v3/include/std/mdspan
>>>>>> @@ -366,6 +366,12 @@ _GLIBCXX_BEGIN_NAMESPACE_VERSION
>>>>>>           class mapping;
>>>>>>       };
>>>>>>
>>>>>> +  struct layout_stride
>>>>>> +  {
>>>>>> +    template<typename _Extents>
>>>>>> +      class mapping;
>>>>>> +  };
>>>>>> +
>>>>>>       namespace __mdspan
>>>>>>       {
>>>>>>         template<typename _Tp>
>>>>>> @@ -434,7 +440,8 @@ _GLIBCXX_BEGIN_NAMESPACE_VERSION
>>>>>>
>>>>>>         template<typename _Mapping>
>>>>>>           concept __standardized_mapping = __mapping_of<layout_left,
>>>> _Mapping>
>>>>>> -                                      || __mapping_of<layout_right,
>>>>>> _Mapping>;
>>>>>> +                                      || __mapping_of<layout_right,
>>>>>> _Mapping>
>>>>>> +                                      || __mapping_of<layout_stride,
>>>>>> _Mapping>;
>>>>>>
>>>>>>         template<typename M>
>>>>>>           concept __mapping_like = requires
>>>>>> @@ -503,6 +510,16 @@ _GLIBCXX_BEGIN_NAMESPACE_VERSION
>>>>>>            : mapping(__other.extents(), __mdspan::__internal_ctor{})
>>>>>>            { }
>>>>>>
>>>>>> +      template<typename _OExtents>
>>>>>> +       requires (is_constructible_v<extents_type, _OExtents>)
>>>>>> +       constexpr explicit(extents_type::rank() > 0)
>>>>>> +       mapping(const layout_stride::mapping<_OExtents>& __other)
>>>>>> +       : mapping(__other.extents(), __mdspan::__internal_ctor{})
>>>>>> +       {
>>>>>> +         __glibcxx_assert(
>>>>>> +           layout_left::mapping<_OExtents>(__other.extents()) ==
>>>> __other);
>>>>>>
>>>>> Could this be *this == other?
>>>>>
>>>>>> +       }
>>>>>> +
>>>>>>           constexpr mapping&
>>>>>>           operator=(const mapping&) noexcept = default;
>>>>>>
>>>>>> @@ -518,8 +535,8 @@ _GLIBCXX_BEGIN_NAMESPACE_VERSION
>>>>>>            constexpr index_type
>>>>>>            operator()(_Indices... __indices) const noexcept
>>>>>>            {
>>>>>> -         return __mdspan::__linear_index_left(
>>>>>> -           this->extents(), static_cast<index_type>(__indices)...);
>>>>>> +         return __mdspan::__linear_index_left(_M_extents,
>>>>>> +           static_cast<index_type>(__indices)...);
>>>>>>            }
>>>>>>
>>>>>>           static constexpr bool
>>>>>> @@ -633,6 +650,16 @@ _GLIBCXX_BEGIN_NAMESPACE_VERSION
>>>>>>            : mapping(__other.extents(), __mdspan::__internal_ctor{})
>>>>>>            { }
>>>>>>
>>>>>> +      template<class _OExtents>
>>>>>> +       requires (is_constructible_v<extents_type, _OExtents>)
>>>>>> +       constexpr explicit(extents_type::rank() > 0)
>>>>>> +       mapping(const layout_stride::mapping<_OExtents>& __other)
>>>> noexcept
>>>>>> +       : mapping(__other.extents(), __mdspan::__internal_ctor{})
>>>>>> +       {
>>>>>> +         __glibcxx_assert(
>>>>>> +           layout_right::mapping<_OExtents>(__other.extents()) ==
>>>>>> __other);
>>>>>>
>>>>> Similary here.
>>>>>
>>>>>> +       }
>>>>>> +
>>>>>>           constexpr mapping&
>>>>>>           operator=(const mapping&) noexcept = default;
>>>>>>
>>>>>> @@ -695,6 +722,192 @@ _GLIBCXX_BEGIN_NAMESPACE_VERSION
>>>>>>            [[no_unique_address]] _Extents _M_extents;
>>>>>>         };
>>>>>>
>>>>>> +  namespace __mdspan
>>>>>> +  {
>>>>>> +    template<typename _Mapping, size_t... _Counts>
>>>>>> +      constexpr typename _Mapping::index_type
>>>>>> +      __offset_impl(const _Mapping& __m, index_sequence<_Counts...>)
>>>>>> noexcept
>>>>>> +      { return __m(((void) _Counts, 0)...); }
>>>>>> +
>>>>>> +    template<typename _Mapping>
>>>>>> +      constexpr typename _Mapping::index_type
>>>>>> +      __offset(const _Mapping& __m) noexcept
>>>>>> +      {
>>>>>>
>>>>> Again, I would define __impl as nested lambda here:
>>>>> auto __impl = [&]<size_t... Countws>(index_seqeunce<_Counts>) noexcept
>>>>> { return  __m(((void) _Counts, 0)...);  }
>>>>>
>>>>>> +       return __offset_impl(__m,
>>>>>> +           make_index_sequence<_Mapping::extents_type::rank()>());
>>>>>> +      }
>>>>>> +
>>>>>> +    template<typename _Mapping, typename... _Indices>
>>>>>> +      constexpr typename _Mapping::index_type
>>>>>> +      __linear_index_strides(const _Mapping& __m,
>>>>>> +                            _Indices... __indices)
>>>>>> +      {
>>>>>> +       using _IndexType = typename _Mapping::index_type;
>>>>>> +       _IndexType __res = 0;
>>>>>> +       if constexpr (sizeof...(__indices) > 0)
>>>>>> +         {
>>>>>> +           auto __update = [&, __pos = 0u](_IndexType __idx) mutable
>>>>>> +             {
>>>>>> +               __res += __idx * __m.stride(__pos++);
>>>>>> +             };
>>>>>> +           (__update(__indices), ...);
>>>>>> +         }
>>>>>> +       return __res;
>>>>>> +      }
>>>>>> +  }
>>>>>> +
>>>>>> +  template<typename _Extents>
>>>>>> +    class layout_stride::mapping
>>>>>> +    {
>>>>>> +      static_assert(__mdspan::__layout_extent<_Extents>,
>>>>>> +       "The size of extents_type must be representable as
>> index_type");
>>>>>> +
>>>>>> +    public:
>>>>>> +      using extents_type = _Extents;
>>>>>> +      using index_type = typename extents_type::index_type;
>>>>>> +      using size_type = typename extents_type::size_type;
>>>>>> +      using rank_type = typename extents_type::rank_type;
>>>>>> +      using layout_type = layout_stride;
>>>>>> +
>>>>>> +      constexpr
>>>>>> +      mapping() noexcept
>>>>>> +      {
>>>>>> +       auto __stride = index_type(1);
>>>>>> +       for(size_t __i = extents_type::rank(); __i > 0; --__i)
>>>>>> +         {
>>>>>> +           _M_strides[__i - 1] = __stride;
>>>>>> +           __stride *= _M_extents.extent(__i - 1);
>>>>>> +         }
>>>>>> +      }
>>>>>> +
>>>>>> +      constexpr
>>>>>> +      mapping(const mapping&) noexcept = default;
>>>>>> +
>>>>>> +      template<__mdspan::__valid_index_type<index_type> _OIndexType>
>>>>>> +       constexpr
>>>>>> +       mapping(const extents_type& __exts,
>>>>>> +               span<_OIndexType, extents_type::rank()> __strides)
>>>> noexcept
>>>>>> +       : _M_extents(__exts)
>>>>>> +       {
>>>>>> +         for(size_t __i = 0; __i < extents_type::rank(); ++__i)
>>>>>> +           _M_strides[__i] = index_type(as_const(__strides[__i]));
>>>>>> +       }
>>>>>> +
>>>>>> +      template<__mdspan::__valid_index_type<index_type> _OIndexType>
>>>>>> +       constexpr
>>>>>> +       mapping(const extents_type& __exts,
>>>>>> +               const array<_OIndexType, extents_type::rank()>&
>>>> __strides)
>>>>>> +       noexcept
>>>>>> +       : mapping(__exts,
>>>>>> +                 span<const _OIndexType,
>>>> extents_type::rank()>(__strides))
>>>>>> +       { }
>>>>>> +
>>>>>> +      template<__mdspan::__mapping_like _StridedMapping>
>>>>>> +       requires (is_constructible_v<extents_type,
>>>>>> +                                    typename
>>>>>> _StridedMapping::extents_type>
>>>>>> +                && _StridedMapping::is_always_unique()
>>>>>> +                && _StridedMapping::is_always_strided())
>>>>>> +       constexpr explicit(!(
>>>>>> +         is_convertible_v<typename _StridedMapping::extents_type,
>>>>>> extents_type>
>>>>>> +         && __mdspan::__standardized_mapping<_StridedMapping>))
>>>>>> +       mapping(const _StridedMapping& __other) noexcept
>>>>>> +       : _M_extents(extents_type(__other.extents()))
>>>>>> +       {
>>>>>> +         if constexpr (extents_type::rank() > 0)
>>>>>> +           for(size_t __i = 0; __i < extents_type::rank(); ++__i)
>>>>>> +             _M_strides[__i] = index_type(__other.stride(__i));
>>>>>> +       }
>>>>>> +
>>>>>> +      constexpr mapping&
>>>>>> +      operator=(const mapping&) noexcept = default;
>>>>>> +
>>>>>> +      constexpr const _Extents&
>>>>>> +      extents() const noexcept { return _M_extents; }
>>>>>> +
>>>>>> +      constexpr array<index_type, extents_type::rank()>
>>>>>> +      strides() const noexcept { return _M_strides; }
>>>>>>
>>>>> Does this constructor work? I mean _M_strides in not std::array,
>>>>> just array. I assumed you would need to create local array and copy
>> them.
>>>>>
>>>>>> +
>>>>>> +      constexpr index_type
>>>>>> +      required_span_size() const noexcept
>>>>>> +      {
>>>>>> +       index_type __ret = 1;
>>>>>> +       for(size_t __i = 0; __i < extents_type::rank(); ++__i)
>>>>>> +         {
>>>>>> +           index_type __ext = _M_extents.extent(__i);
>>>>>> +           if(__ext == 0)
>>>>>> +             return 0;
>>>>>> +           __ret += (__ext - 1) * _M_strides[__i];
>>>>>> +         }
>>>>>> +       return __ret;
>>>>>> +      }
>>>>>> +
>>>>>> +      template<__mdspan::__valid_index_type<index_type>... _Indices>
>>>>>> +       requires (sizeof...(_Indices) == extents_type::rank())
>>>>>> +       constexpr index_type operator()(_Indices... __indices) const
>>>>>> noexcept
>>>>>> +       {
>>>>>> +         return __mdspan::__linear_index_strides(*this,
>>>>>> +           static_cast<index_type>(__indices)...);
>>>>>> +       }
>>>>>> +
>>>>>> +      static constexpr bool
>>>>>> +      is_always_unique() noexcept { return true; }
>>>>>> +
>>>>>> +      static constexpr bool
>>>>>> +      is_always_exhaustive() noexcept { return false; }
>>>>>> +
>>>>>> +      static constexpr bool
>>>>>> +      is_always_strided() noexcept { return true; }
>>>>>> +
>>>>>> +      static constexpr bool
>>>>>> +      is_unique() noexcept { return true; }
>>>>>> +
>>>>>> +      constexpr bool
>>>>>> +      is_exhaustive() const noexcept
>>>>>> +      {
>>>>>> +       if constexpr (extents_type::rank() == 0)
>>>>>> +         return true;
>>>>>> +       else
>>>>>> +       {
>>>>>> +         // Under the assumption that the mapping is unique and has
>>>>>> positive
>>>>>> +         // size, the mapping is exhaustive, if and only if the
>> largest
>>>>>> value
>>>>>> +         // returned by m(i...) is within the required range.
>>>>>> +         // However, the standard requires implementing a condition
>>>>>> that's not
>>>>>> +         // always true when the size of the mapping is zero.
>>>>>>
>>>>> With the precondition on the constructor
>>>>> https://eel.is/c++draft/mdspan.layout.stride#cons-4.3,
>>>>> I do not think you can produce a layout_stride where this would
>> produce a
>>>>> different result.
>>>>> If you agree I would change the comment that "C++23
>>>>> [mdspan.layout.stride.cons] p4.3"
>>>>> implies that this produces the same result.
>>>>
>>>> If the size of the extent is zero, then almost any choice of strides
>>>> will work as a counter example. For example:
>>>>
>>>>      extents = [0, 1, 1]
>>>>      strides = [2, 4, 8]
>>>>
>>>> is exhaustive, because every element (of the empty set) is covered.
>>>> However, [2, 4, 8] don't match the formula in the spec [1].
>>>>
>>> My point is that constructing layout_stride with extents and strides from
>>> the
>>> example violates the precondition on the constructor:
>>> https://eel.is/c++draft/mdspan.layout.stride#cons-4.3.
>>> So this function cannot observe the difference.
>>>
>>
>> I don't see it. I'll run through it step by step.
>>
>> Pick `p[i] = i` as the permutation and check that:
>>      i == 1: s[p[1]] >= s[p[0]] * e[p[0]]
>>                 s[1] >= s[0] * e[0]
>>                   4  >= 2 * 0
>>
>>      i == 2:    s[2] >= s[1] * e[1]
>>                   8  >= 4 * 1
>>
>> Therefore, rank_ == 3 > 0 and there exists a permutation P
>> such that
>>       s[P[i]] >= s[P[i-1]] * e[P[i-1]
>> for all i in [1, rank_).
>>
>> I think that covers the entire precondition.
>>
> Indeed, however if we consider following:
>    extents = [1, 3, 0]
>    strides = [2, 8, 8]
> This is exhaustive (and your checks answer that), and technically
> is_strided.
> However, the constructor will not accept it. I think this is general
> question on
> empty mappings and layout strided.
> 

Thank you, I see what you mean even if this precise case would
work with P = [0, 2, 1]:

   i == 1: s[p[1]] >= s[p[0]] * e[p[0]]
                 8 >= 2 * 1

   i == 2: s[p[2]] >= s[p[1]] * e[p[1]]
                 8 >= 0 * 8

Which leads us to:

[Note 1: For layout_stride, this condition is necessary and sufficient
  for is_unique() to be true. — end note]

The "nesessary" is, I believe, incorrect, because of this variant of
the counter example you provide:

   extents = [2, 2, 2, 0]
   strides = [2, 2, 2, 2]

> 
>>>>
>>>> If the size of extents isn't zero, then the two approaches are the
>>>> same.
>>>>
>>>> [1]: https://eel.is/c++draft/mdspan.layout.stride#obs-5.2
>>>>
>>>>>>
>>>>>> +         auto __size = __mdspan::__fwd_prod(_M_extents,
>>>>>> extents_type::rank());
>>>>>> +         if (__size == 0)
>>>>>> +           return true;
>>>>>> +         return __size == required_span_size();
>>>>>> +       }
>>>>>> +      }
>>>>>> +
>>>>>> +      static constexpr bool
>>>>>> +      is_strided() noexcept { return true; }
>>>>>> +
>>>>>> +      constexpr index_type
>>>>>> +      stride(rank_type __r) const noexcept { return _M_strides[__r];
>> }
>>>>>> +
>>>>>> +      template<__mdspan::__mapping_like _OMapping>
>>>>>> +       requires ((extents_type::rank() ==
>>>> _OMapping::extents_type::rank())
>>>>>> +                 && _OMapping::is_always_strided())
>>>>>> +       friend constexpr bool
>>>>>> +       operator==(const mapping& __self, const _OMapping& __other)
>>>>>> noexcept
>>>>>> +       {
>>>>>> +         if(__self.extents() != __other.extents())
>>>>>> +           return false;
>>>>>> +         if constexpr (extents_type::rank() > 0)
>>>>>> +           for(size_t __i = 0; __i < extents_type::rank(); ++__i)
>>>>>> +             if (__self.stride(__i) != __other.stride(__i))
>>>>>> +               return false;
>>>>>> +         return __mdspan::__offset(__other) == 0;
>>>>>> +       }
>>>>>> +
>>>>>> +    private:
>>>>>> +      using _S_strides_t = typename __array_traits<index_type,
>>>>>> +
>>>>>> extents_type::rank()>::_Type;
>>>>>> +      [[no_unique_address]] _Extents _M_extents;
>>>>>> +      [[no_unique_address]] _S_strides_t _M_strides;
>>>>>> +    };
>>>>>> +
>>>>>>     _GLIBCXX_END_NAMESPACE_VERSION
>>>>>>     }
>>>>>>     #endif
>>>>>> --
>>>>>> 2.49.0
>>>>>>
>>>>>>
>>>>>
>>>>
>>>>
>>>
>>
>>
> 



More information about the Libstdc++ mailing list