[PATCH v1 3/6] libstdc++: Improve low-rank layout_{left, right}::stride.

Luc Grosheintz luc.grosheintz@gmail.com
Wed Jul 30 08:56:29 GMT 2025



On 7/28/25 13:04, Tomasz Kaminski wrote:
> On Mon, Jul 28, 2025 at 10:24 AM Tomasz Kaminski <tkaminsk@redhat.com>
> wrote:
> 
>>
>>
>>
>> On Mon, Jul 28, 2025 at 10:03 AM Luc Grosheintz <luc.grosheintz@gmail.com>
>> wrote:
>>
>>>
>>>
>>> On 7/28/25 08:02, Tomasz Kaminski wrote:
>>>> On Sun, Jul 27, 2025 at 2:47 PM Luc Grosheintz <
>>> luc.grosheintz@gmail.com>
>>>> wrote:
>>>>
>>>>> The methods layout_{left,right}::mapping::stride are defined
>>>>> as
>>>>>
>>>>>     \prod_{i = 0}^r E[i]
>>>>>     \prod_{i = r+1}^n E[i]
>>>>>
>>>>> This is computed as the product of a pre-comupted static product and
>>> the
>>>>> product of the required dynamic extents.
>>>>>
>>>>> Disassembly shows that even for low-rank extents, i.e. rank == 1 and
>>>>> rank == 2, with at least one  dynamic extent, the generated code loads
>>>>> two values; and then runs the loop over at most one element, e.g.
>>>>>
>>>>>    220:  48 8b 0c f5 00 00 00   mov    rcx,QWORD PTR [rsi*8+0x0]
>>>>>    227:  00
>>>>>    228:  48 8b 04 f5 00 00 00   mov    rax,QWORD PTR [rsi*8+0x0]
>>>>>    22f:  00
>>>>>    230:  48 c1 e1 02            shl    rcx,0x2
>>>>>    234:  74 1a                  je     250 <stride_left_d5+0x30>
>>>>>    236:  48 01 f9               add    rcx,rdi
>>>>>    239:  0f 1f 80 00 00 00 00   nop    DWORD PTR [rax+0x0]
>>>>>    240:  48 63 17               movsxd rdx,DWORD PTR [rdi]
>>>>>    243:  48 83 c7 04            add    rdi,0x4
>>>>>    247:  48 0f af c2            imul   rax,rdx
>>>>>    24b:  48 39 f9               cmp    rcx,rdi
>>>>>    24e:  75 f0                  jne    240 <stride_left_d5+0x20>
>>>>>    250:  c3                     ret
>>>>>
>>>>> If there's no dynamic extents, it simply loads the precomputed product
>>>>> of static extents.
>>>>>
>>>>> For rank == 1 the answer is constant `1`; for rank == 2 it's either 1
>>> or
>>>>> extents.extent(k), with k == 0 for layout_left and k == 1 for
>>>>> layout_right.
>>>>>
>>>>> Consider,
>>>>>
>>>>>     using Ed = std::extents<int, dyn>;
>>>>>     int stride_left_d(const std::layout_left::mapping<Ed>& m, size_t r)
>>>>>     { return m.stride(r); }
>>>>>
>>>>>     using E3d = std::extents<int, 3, dyn>;
>>>>>     int stride_left_3d(const std::layout_left::mapping<E3d>& m, size_t
>>> r)
>>>>>     { return m.stride(r); }
>>>>>
>>>>>     using Ed5 = std::extents<int, dyn, 5>;
>>>>>     int stride_left_d5(const std::layout_left::mapping<Ed5>& m, size_t
>>> r)
>>>>>     { return m.stride(r); }
>>>>>
>>>>> The optimized code for these three cases is:
>>>>>
>>>>>     0000000000000060 <stride_left_d>:
>>>>>     60:  b8 01 00 00 00         mov    eax,0x1
>>>>>     65:  c3                     ret
>>>>>
>>>>>     0000000000000090 <stride_left_3d>:
>>>>>     90:  48 83 fe 01            cmp    rsi,0x1
>>>>>     94:  19 c0                  sbb    eax,eax
>>>>>     96:  83 e0 fe               and    eax,0xfffffffe
>>>>>     99:  83 c0 03               add    eax,0x3
>>>>>     9c:  c3                     ret
>>>>>
>>>>>     00000000000000a0 <stride_left_d5>:
>>>>>     a0:  b8 01 00 00 00         mov    eax,0x1
>>>>>     a5:  48 85 f6               test   rsi,rsi
>>>>>     a8:  74 02                  je     ac <stride_left_d5+0xc>
>>>>>     aa:  8b 07                  mov    eax,DWORD PTR [rdi]
>>>>>    i ac:  c3                     ret
>>>>>
>>>>> For rank == 1 it simply returns 1 (as expected). For rank == 2, it
>>>>> either implements a branchless formula, or conditionally loads one
>>>>> value. In all cases involving a dynamic extent this seems like it's
>>>>> always doing clearly less work, both in terms of computation and loads.
>>>>>
>>>>> For rank == 2, it trades loading one value for a branchless sequence of
>>>>> four instructions that don't require loading any values.
>>>>>
>>>> I will put this optimization into the __fwd_prod and __rev_pord
>>> functions,
>>>> so it will be applied for all uses. This will also avoid us creating
>>> this
>>>> caching
>>>> tables for such small ranks.
>>>
>>> The problem is that we don't have the same amount of information in
>>> the stride and __fwd_prod. The valid values for __r are 1, ..., rank;
>>> for __fwd_prod it's inclusive, while in stride it's exclusive. Therefore,
>>> we can do the optimization with one comparison less in stride than in
>>> __fwd_prod.
>>>
>> Make sense, and I am OK having this optimization there. However, out of
>> curiosity, if we put always_inline on the __fwd_prod, wouldn't compiler be
>> able
>> to eliminate __rank == __i. Also once we have separate call to __size, we
>> could
>> put assertions (just as a comment) that for __fwd_prod and __rev_prod __r
>> is
>> required to be smaller than rank().
>>
> I think I would pursue that direction, and have similar if constexpr in the
> __size function:
> __rank == 0 -> 1
> __rank == 1 -> extent(0)
> __rank == 2 -> extent(0) * extent(1)
> __rank == rank_dynamic() -> __ext_prod with __sta_ext == 1
> otherwise -> use __sta_extr.

We know that __size is optimized very aggressively and effectively. On
-O2 everything is inlined and no symbols are created, no static storage
is needed to compute __size. This is regardless of rank.

The following:

   using E1 = std::extents<int, dyn>;
   int size1(const std::mdspan<double, E1>& md)
   { return md.size(); }

   using E2 = std::extents<int, 3, 5, dyn, dyn, 7, dyn, 11>;
   int size2(const std::mdspan<double, E2>& md)
   { return md.size(); }

   using E3 = std::extents<int, 3, dyn>;
   int size3(const std::mdspan<double, E3>& md)
   { return md.size(); }

when compiled with -O2 -c size.o

$ nm size.o
0000000000000000 T size1
0000000000000010 T size2
0000000000000030 T size3

$ objdump -h size.o
size.o:     file format elf64-x86-64

Sections:
Idx Name          Size      VMA               LMA               File off  Algn
   0 .text         00000037  0000000000000000  0000000000000000  00000040  2**4
                   CONTENTS, ALLOC, LOAD, READONLY, CODE
   1 .data         00000000  0000000000000000  0000000000000000  00000077  2**0
                   CONTENTS, ALLOC, LOAD, DATA
   2 .bss          00000000  0000000000000000  0000000000000000  00000077  2**0
                   ALLOC
   3 .comment      0000002b  0000000000000000  0000000000000000  00000077  2**0
                   CONTENTS, READONLY
   4 .note.GNU-stack 00000000  0000000000000000  0000000000000000  000000a2  2**0
                   CONTENTS, READONLY
   5 .note.gnu.property 00000030  0000000000000000  0000000000000000  000000a8  2**3
                   CONTENTS, ALLOC, LOAD, READONLY, DATA
   6 .eh_frame     00000058  0000000000000000  0000000000000000  000000d8  2**3
                   CONTENTS, ALLOC, LOAD, RELOC, READONLY, DATA

$ objdump -d size.o
0000000000000000 <size1>:
    0:	8b 07                	mov    (%rdi),%eax
    2:	c3                   	ret
    3:	66 66 2e 0f 1f 84 00 	data16 cs nopw 0x0(%rax,%rax,1)
    a:	00 00 00 00
    e:	66 90                	xchg   %ax,%ax

0000000000000010 <size2>:
   10:	48 63 07             	movslq (%rdi),%rax
   13:	48 63 57 04          	movslq 0x4(%rdi),%rdx
   17:	48 0f af c2          	imul   %rdx,%rax
   1b:	0f af 47 08          	imul   0x8(%rdi),%eax
   1f:	69 c0 83 04 00 00    	imul   $0x483,%eax,%eax
   25:	c3                   	ret
   26:	66 2e 0f 1f 84 00 00 	cs nopw 0x0(%rax,%rax,1)
   2d:	00 00 00

0000000000000030 <size3>:
   30:	48 63 07             	movslq (%rdi),%rax
   33:	8d 04 40             	lea    (%rax,%rax,2),%eax
   36:	c3                   	ret

Hence, I'm not sure where I need to look to find the cost of not putting
the constexpr-if in __size. If I can see why it's slower, I'm happy to
add them. However, since I've looked and can't see a cost, I'm reluctant
to add them because they require test coverage for each branch (otherwise
they're never compiled and could contain complete non-sense).

> 
>>
>>> I purposefully checked mdspan::size (part of the previous commit); and
>>> on optimized builds it fully unrolls the loop and does everything
>>> automatically, meaning it doesn't need help and we're not repeating
>>> the same optimization several times.
>>>
>>> As for avoiding the tables for small ranks, we can refactor _RevProd
>>> as follows:
>>>
>>>     template<size_t... _Extents>
>>>       struct _RevProd
>>>       {
>>>          size_t _S_value(size_t i)
>>>          { return _S_data[i]; }
>>>
>>>       private:
>>>          array _S_data = consteval ...;
>>>       }
>>>
>>> and use partial specialization to eliminate _S_data.
>>>
>>>>
>>>>>
>>>>> libstdc++-v3/ChangeLog:
>>>>>
>>>>>           * include/std/mdspan (layout_left::mapping::stride): Optimize
>>>>>           for rank <= 2.
>>>>>           (layout_right::mapping::stride): Ditto.
>>>>>
>>>>> Signed-off-by: Luc Grosheintz <luc.grosheintz@gmail.com>
>>>>> ---
>>>>>    libstdc++-v3/include/std/mdspan | 14 ++++++++++++--
>>>>>    1 file changed, 12 insertions(+), 2 deletions(-)
>>>>>
>>>>> diff --git a/libstdc++-v3/include/std/mdspan
>>>>> b/libstdc++-v3/include/std/mdspan
>>>>> index 06ccf3e3827..f288af96cdb 100644
>>>>> --- a/libstdc++-v3/include/std/mdspan
>>>>> +++ b/libstdc++-v3/include/std/mdspan
>>>>> @@ -652,7 +652,12 @@ _GLIBCXX_BEGIN_NAMESPACE_VERSION
>>>>>          requires (extents_type::rank() > 0)
>>>>>          {
>>>>>           __glibcxx_assert(__i < extents_type::rank());
>>>>> -       return __mdspan::__fwd_prod(_M_extents, __i);
>>>>> +       if constexpr (extents_type::rank() == 1)
>>>>> +         return 1;
>>>>> +       else if constexpr (extents_type::rank() == 2)
>>>>> +         return __i == 0 ? 1 : _M_extents.extent(0);
>>>>> +       else
>>>>> +         return __mdspan::__fwd_prod(_M_extents, __i);
>>>>>          }
>>>>>
>>>>>          template<typename _OExtents>
>>>>> @@ -797,7 +802,12 @@ _GLIBCXX_BEGIN_NAMESPACE_VERSION
>>>>>          requires (extents_type::rank() > 0)
>>>>>          {
>>>>>           __glibcxx_assert(__i < extents_type::rank());
>>>>> -       return __mdspan::__rev_prod(_M_extents, __i);
>>>>> +       if constexpr (extents_type::rank() == 1)
>>>>> +         return 1;
>>>>> +       else if constexpr (extents_type::rank() == 2)
>>>>> +         return __i == 0 ? _M_extents.extent(1) : 1;
>>>>> +       else
>>>>> +         return __mdspan::__rev_prod(_M_extents, __i);
>>>>>          }
>>>>>
>>>>>          template<typename _OExtents>
>>>>> --
>>>>> 2.50.0
>>>>>
>>>>>
>>>>
>>>
>>>
> 



More information about the Libstdc++ mailing list