Vir Campestris wrote:

> It is a microbenchmark.
> I was trying to understand cache performance on my system.
> Over in comp.lang.c++ you'll find a thread about Sieve or Eratosthenes.
>
> I and another poster are trying to optimise it. Not for any good reason
> of course... but I was just curious about some of the results I was
> getting.

A few years ago this guy popped-up on youtube, saying "I was the author
of Task Manager on Windows NT", I ignored the recommendation for quite a
while thinking he just wanted people to blow smoke up his arse,
eventually I watched and he seems a decent guy ... anyway he did a drag
racing series using sieve prime algorithm

<https://www.youtube.com/playlist?list=PLF2KJ6Gy3cZ5Er-1eF9fN1Hgw_xkoD9V1>

Others have contributed in all manner of languages

<https://github.com/PlummersSoftwareLLC/Primes>

Vir Campestris wrote:

> The code I have put in a signature block; there's no point in risking
> someone posting it again. I've commented it, but no doubt not in all the
> right places! I'd be interested to know what results other people get.

Built and run under Win11 with g++ 13.2.0
on a laptop with
i7-11370H CPU (3.3GHz base, 4.28GHz turbo)
L1 d-cache 48kB per core
L2 cache 1.2MB
L3 cache 12MB
16GB memory

Size 1 gave 5.12047e+09 bytes/second.
Size 2 gave 5.76467e+09 bytes/second.
Size 4 gave 5.78979e+09 bytes/second.
Size 8 gave 5.6642e+09 bytes/second.
Size 16 gave 5.55303e+09 bytes/second.
Size 32 gave 4.4151e+09 bytes/second.
Size 64 gave 4.67783e+09 bytes/second.
Size 128 gave 4.85115e+09 bytes/second.
Size 256 gave 4.88147e+09 bytes/second.
Size 512 gave 4.5214e+09 bytes/second.
Size 1024 gave 4.64794e+09 bytes/second.
Size 2048 gave 4.69149e+09 bytes/second.
Size 4096 gave 4.69416e+09 bytes/second.
Size 8192 gave 4.73425e+09 bytes/second.
Size 16384 gave 4.76353e+09 bytes/second.
Size 32768 gave 4.70864e+09 bytes/second.
Size 65536 gave 4.75244e+09 bytes/second.
Size 131072 gave 4.76452e+09 bytes/second.
Size 262144 gave 4.69839e+09 bytes/second.
Size 524288 gave 4.71094e+09 bytes/second.
Size 1048576 gave 4.62922e+09 bytes/second.
Size 2097152 gave 4.55435e+09 bytes/second.
Size 4194304 gave 4.58315e+09 bytes/second.
Size 8388608 gave 4.55917e+09 bytes/second.
Size 16777216 gave 4.60661e+09 bytes/second.
Size 33554432 gave 4.63509e+09 bytes/second.
Size 67108864 gave 4.62349e+09 bytes/second.
Size 134217728 gave 4.62406e+09 bytes/second.

Andy Burns wrote:

> Size 1 gave 5.12047e+09 bytes/second.
> Size 2 gave 5.76467e+09 bytes/second.
> Size 4 gave 5.78979e+09 bytes/second.
> Size 8 gave 5.6642e+09 bytes/second.
> Size 16 gave 5.55303e+09 bytes/second.
> Size 32 gave 4.4151e+09 bytes/second.
> Size 64 gave 4.67783e+09 bytes/second.
> Size 128 gave 4.85115e+09 bytes/second.
> Size 256 gave 4.88147e+09 bytes/second.
> Size 512 gave 4.5214e+09 bytes/second.
> Size 1024 gave 4.64794e+09 bytes/second.
> Size 2048 gave 4.69149e+09 bytes/second.
> Size 4096 gave 4.69416e+09 bytes/second.
> Size 8192 gave 4.73425e+09 bytes/second.
> Size 16384 gave 4.76353e+09 bytes/second.
> Size 32768 gave 4.70864e+09 bytes/second.
> Size 65536 gave 4.75244e+09 bytes/second.
> Size 131072 gave 4.76452e+09 bytes/second.
> Size 262144 gave 4.69839e+09 bytes/second.
> Size 524288 gave 4.71094e+09 bytes/second.
> Size 1048576 gave 4.62922e+09 bytes/second.
> Size 2097152 gave 4.55435e+09 bytes/second.
> Size 4194304 gave 4.58315e+09 bytes/second.
> Size 8388608 gave 4.55917e+09 bytes/second.
> Size 16777216 gave 4.60661e+09 bytes/second.
> Size 33554432 gave 4.63509e+09 bytes/second.
> Size 67108864 gave 4.62349e+09 bytes/second.
> Size 134217728 gave 4.62406e+09 bytes/second.

and (perversely?) with -Ofast

Size 1 gave 4.64966e+09 bytes/second.
Size 2 gave 9.39481e+09 bytes/second.
Size 4 gave 1.59656e+10 bytes/second.
Size 8 gave 2.6789e+10 bytes/second.
Size 16 gave 4.42999e+10 bytes/second.
Size 32 gave 4.51076e+10 bytes/second.
Size 64 gave 4.90103e+10 bytes/second.
Size 128 gave 4.08414e+10 bytes/second.
Size 256 gave 4.36978e+10 bytes/second.
Size 512 gave 4.54119e+10 bytes/second.
Size 1024 gave 4.50594e+10 bytes/second.
Size 2048 gave 4.71622e+10 bytes/second.
Size 4096 gave 4.60261e+10 bytes/second.
Size 8192 gave 3.84764e+10 bytes/second.
Size 16384 gave 3.90878e+10 bytes/second.
Size 32768 gave 3.76718e+10 bytes/second.
Size 65536 gave 3.87339e+10 bytes/second.
Size 131072 gave 3.89058e+10 bytes/second.
Size 262144 gave 3.85662e+10 bytes/second.
Size 524288 gave 4.19209e+10 bytes/second.
Size 1048576 gave 3.00962e+10 bytes/second.
Size 2097152 gave 1.37616e+10 bytes/second.
Size 4194304 gave 1.4136e+10 bytes/second.
Size 8388608 gave 1.39413e+10 bytes/second.
Size 16777216 gave 1.35747e+10 bytes/second.
Size 33554432 gave 1.37286e+10 bytes/second.
Size 67108864 gave 1.38163e+10 bytes/second.
Size 134217728 gave 1.29987e+10 bytes/second.

Andy Burns wrote:

> Vir Campestris wrote:

>> The code I have put in a signature block; there's no point in risking
>> someone posting it again. I've commented it, but no doubt not in all the
>> right places! I'd be interested to know what results other people get.

> Built and run under Win11 with g++ 13.2.0
> on a laptop with
> i7-11370H CPU (3.3GHz base, 4.28GHz turbo)
> L1 d-cache 48kB per core
> L2 cache 1.2MB
> L3 cache 12MB
> 16GB memory

> Size 1 gave 5.12047e+09 bytes/second.
> Size 2 gave 5.76467e+09 bytes/second.
> Size 4 gave 5.78979e+09 bytes/second.
> Size 8 gave 5.6642e+09 bytes/second.
> Size 16 gave 5.55303e+09 bytes/second.
> Size 32 gave 4.4151e+09 bytes/second.
> Size 64 gave 4.67783e+09 bytes/second.
> Size 128 gave 4.85115e+09 bytes/second.
> Size 256 gave 4.88147e+09 bytes/second.
> Size 512 gave 4.5214e+09 bytes/second.
> Size 1024 gave 4.64794e+09 bytes/second.
> Size 2048 gave 4.69149e+09 bytes/second.
> Size 4096 gave 4.69416e+09 bytes/second.
> Size 8192 gave 4.73425e+09 bytes/second.
> Size 16384 gave 4.76353e+09 bytes/second.
> Size 32768 gave 4.70864e+09 bytes/second.
> Size 65536 gave 4.75244e+09 bytes/second.
> Size 131072 gave 4.76452e+09 bytes/second.
> Size 262144 gave 4.69839e+09 bytes/second.
> Size 524288 gave 4.71094e+09 bytes/second.
> Size 1048576 gave 4.62922e+09 bytes/second.
> Size 2097152 gave 4.55435e+09 bytes/second.
> Size 4194304 gave 4.58315e+09 bytes/second.
> Size 8388608 gave 4.55917e+09 bytes/second.
> Size 16777216 gave 4.60661e+09 bytes/second.
> Size 33554432 gave 4.63509e+09 bytes/second.
> Size 67108864 gave 4.62349e+09 bytes/second.
> Size 134217728 gave 4.62406e+09 bytes/second.

This data set has a 5/4 ratio best to worst with noise.

Andy Burns wrote:

> and (perversely?) with -Ofast

I must look at the exponent, not just the mantissa
I must look at the exponent, not just the mantissa
I must look at the exponent, not just the mantissa

Andy Burns wrote:

> Andy Burns wrote:

>> Size 1 gave 5.12047e+09 bytes/second.
>> Size 2 gave 5.76467e+09 bytes/second.
>> Size 4 gave 5.78979e+09 bytes/second.
>> Size 8 gave 5.6642e+09 bytes/second.
>> Size 16 gave 5.55303e+09 bytes/second.
>> Size 32 gave 4.4151e+09 bytes/second.
>> Size 64 gave 4.67783e+09 bytes/second.
>> Size 128 gave 4.85115e+09 bytes/second.
>> Size 256 gave 4.88147e+09 bytes/second.
>> Size 512 gave 4.5214e+09 bytes/second.
>> Size 1024 gave 4.64794e+09 bytes/second.
>> Size 2048 gave 4.69149e+09 bytes/second.
>> Size 4096 gave 4.69416e+09 bytes/second.
>> Size 8192 gave 4.73425e+09 bytes/second.
>> Size 16384 gave 4.76353e+09 bytes/second.
>> Size 32768 gave 4.70864e+09 bytes/second.
>> Size 65536 gave 4.75244e+09 bytes/second.
>> Size 131072 gave 4.76452e+09 bytes/second.
>> Size 262144 gave 4.69839e+09 bytes/second.
>> Size 524288 gave 4.71094e+09 bytes/second.
>> Size 1048576 gave 4.62922e+09 bytes/second.
>> Size 2097152 gave 4.55435e+09 bytes/second.
>> Size 4194304 gave 4.58315e+09 bytes/second.
>> Size 8388608 gave 4.55917e+09 bytes/second.
>> Size 16777216 gave 4.60661e+09 bytes/second.
>> Size 33554432 gave 4.63509e+09 bytes/second.
>> Size 67108864 gave 4.62349e+09 bytes/second.
>> Size 134217728 gave 4.62406e+09 bytes/second.

> and (perversely?) with -Ofast
{ A
> Size 1 gave 4.64966e+09 bytes/second.
> Size 2 gave 9.39481e+09 bytes/second.
> Size 4 gave 1.59656e+10 bytes/second.
> Size 8 gave 2.6789e+10 bytes/second.
> Size 16 gave 4.42999e+10 bytes/second.
} A
{ B
> Size 32 gave 4.51076e+10 bytes/second.
> Size 64 gave 4.90103e+10 bytes/second.
> Size 128 gave 4.08414e+10 bytes/second.
> Size 256 gave 4.36978e+10 bytes/second.
> Size 512 gave 4.54119e+10 bytes/second.
> Size 1024 gave 4.50594e+10 bytes/second.
> Size 2048 gave 4.71622e+10 bytes/second.
> Size 4096 gave 4.60261e+10 bytes/second.
} B
{ C
> Size 8192 gave 3.84764e+10 bytes/second.
> Size 16384 gave 3.90878e+10 bytes/second.
> Size 32768 gave 3.76718e+10 bytes/second.
> Size 65536 gave 3.87339e+10 bytes/second.
> Size 131072 gave 3.89058e+10 bytes/second.
> Size 262144 gave 3.85662e+10 bytes/second.
> Size 524288 gave 4.19209e+10 bytes/second.
} C
{ D
> Size 1048576 gave 3.00962e+10 bytes/second.
> Size 2097152 gave 1.37616e+10 bytes/second.
> Size 4194304 gave 1.4136e+10 bytes/second.
> Size 8388608 gave 1.39413e+10 bytes/second.
> Size 16777216 gave 1.35747e+10 bytes/second.
> Size 33554432 gave 1.37286e+10 bytes/second.
> Size 67108864 gave 1.38163e+10 bytes/second.
> Size 134217728 gave 1.29987e+10 bytes/second.
} D

A displays BW follows access size

B displays essentially the same performance throughout

C displays a minor step down in performance due to L2 or TLB effects 5 to 4 ratio

D displays a major stepdown in performance at DIV 3 compared to C 3 to 1 ratio

Andy Burns wrote:
> Vir Campestris wrote:
>
>> It is a microbenchmark.
>> I was trying to understand cache performance on my system.
>> Over in comp.lang.c++ you'll find a thread about Sieve or Eratosthenes.
>>
>> I and another poster are trying to optimise it. Not for any good
>> reason of course... but I was just curious about some of the results I
>> was getting.
>
> A few years ago this guy popped-up on youtube, saying "I was the author
> of Task Manager on Windows NT", I ignored the recommendation for quite a
> while thinking he just wanted people to blow smoke up his arse,
> eventually I watched and he seems a decent guy ... anyway he did a drag
> racing series using sieve prime algorithm
>
> <https://www.youtube.com/playlist?list=PLF2KJ6Gy3cZ5Er-1eF9fN1Hgw_xkoD9V1>

Dave is a Good Guy!

I have contributed to his "smallest windows program ever" with a small
tweak that saved 4 (or 6?) bytes.

The sieve algorithm used above is more advanced than the fixed modulo 30
program I wrote a ffew decades ago, but similar in approach.

I decided to use mod 30 because as soon as you get to 30+, there are
exactly 8 possible prime locations in each block, so it fits perfectly
in a byte map. :-)

1,7,11,13,17,19,23,29
>
> Others have contributed in all manner of languages
>
> <https://github.com/PlummersSoftwareLLC/Primes>

Terje

--
- <Terje.Mathisen at tmsw.no>
"almost all programming can be viewed as an exercise in caching"

On 12/30/2023 7:57 AM, Terje Mathisen wrote:
> Andy Burns wrote:
>> Vir Campestris wrote:
>>
>>> It is a microbenchmark.
>>> I was trying to understand cache performance on my system.
>>> Over in comp.lang.c++ you'll find a thread about Sieve or Eratosthenes.
>>>
>>> I and another poster are trying to optimise it. Not for any good
>>> reason of course... but I was just curious about some of the results
>>> I was getting.
>>
>> A few years ago this guy popped-up on youtube, saying "I was the
>> author of Task Manager on Windows NT", I ignored the recommendation
>> for quite a while thinking he just wanted people to blow smoke up his
>> arse, eventually I watched and he seems a decent guy ... anyway he did
>> a drag racing series using sieve prime algorithm
>>
>> <https://www.youtube.com/playlist?list=PLF2KJ6Gy3cZ5Er-1eF9fN1Hgw_xkoD9V1>
>
> Dave is a Good Guy!
>

Yes, and generally interesting videos on YouTube.

He did interviews with David Cutler and Raymond Chen as well, both
fairly big names at Microsoft, which were also kinda interesting.

> I have contributed to his "smallest windows program ever" with a small
> tweak that saved 4 (or 6?) bytes.
>
> The sieve algorithm used above is more advanced than the fixed modulo 30
> program I wrote a ffew decades ago, but similar in approach.
>
> I decided to use mod 30 because as soon as you get to 30+, there are
> exactly 8 possible prime locations in each block, so it fits perfectly
> in a byte map. :-)
>
> 1,7,11,13,17,19,23,29

I before used a prime sieve to test static vs dynamic types performance
in my compiler, and noted that the relative performance difference was
smaller than might otherwise be expected (intuitively, one would expect
dynamic types to be horridly slow, since nearly every operation is
effectively a runtime call).

Namely, I only saw around a 3x difference at the time, rather than
around a 10x+ difference, which is more what I would expect.

Granted, a prime sieve is also not exactly something that performs well
on my ISA design (but, slightly less "aggressively badly" than a the
"recursive Fibonacci number" algorithm, *).

*, IIRC:
long rfib(long x)
{ return((x<2)?1:(rfib(x-1)+rfib(x-2))); }
Which is basically a stress test of function call and return performance...

>>
>> Others have contributed in all manner of languages
>>
>> <https://github.com/PlummersSoftwareLLC/Primes>
>
>
> Terje
>
>

BGB wrote:

> long rfib(long x)
> { return((x<2)?1:(rfib(x-1)+rfib(x-2))); }

At first I thought the following looked pretty good:

rfib:
ENTER R29,R0,#0 // just 3 preserved registers.

CMP R2,R1,#2
PGE R2,TT
MOV R1,#1 // return 1
EXIT R29,R0,#0

SUB R30,R1,#2
SUB R1,R1,#1
CALL rfib // rfib(n-1)
MOV R29,R1
MOV R1,R30
CALL rfib // rfib(n-2)

ADD R1,R1,R29
EXIT R29,R0,#0 // return rfib(n-1)+rfib(n-2)

But then I moved save and restore to the rfib() calls making the last rung
on the call tree less expensive by 6 memory references each, and by saving
restoring only 2 registers rather than 3 at the cost of a MOV, so each non-
leaf level in the call tree saves/restores only 4 registers instead of 6::

rfib:
SUB R2,R1,#2
PNLT R2,TT // (n-2) < 0
MOV R1,#1 // return 1
RET

ENTER R30,R0,#0 // just 2 preserved registers.
MOV R30,R2
SUB R1,R1,#1
CALL rfib // rfib(n-1)
MOV R2,R1
MOV R1,R30
MOV R30,R2
CALL rfib // rfib(n-2)

ADD R1,R1,R30
EXIT R30,R0,#0 // return rfib(n-1)+rfib(n-2)

But then I realized that rfib(n-2) has already been computed by rfib(n-1).
Since the size of the redundant element on the stack is constant, rfib(n-1)
makes a location available to rfib(n-2) so only 1 call is required instead
of 2::

rfib:
ENTER R0,R0,#8
STD #0,[SP+8] // setup rfib[n-2]
CALL rfib1 // recurse through rfib1
EXIT R0,R0,#8

rfib1:
SUB R2,R1,#2
PNLT R2,TT // (n-2) < 0
MOV R1,#1 // return 1
RET

ENTER R0,R0,#8 // no preserved registers just return address

SUB R1,R1,#1
CALL rfib1 // rfib(n-1)
LDD R2,[SP] // save rfib(n-2)

ADD R1,R1,R30
STD R1,[SP+16] // restore rfib(n-2)
EXIT R0,R0,#8 // return rfib(n-1)+rfib(n-2)

.....

On 12/30/2023 4:05 PM, MitchAlsup wrote:
> BGB wrote:
>
>>    long rfib(long x)
>>      { return((x<2)?1:(rfib(x-1)+rfib(x-2))); }
>
> At first I thought the following looked pretty good:
>
> rfib:
>      ENTER   R29,R0,#0          // just 3 preserved registers.
>
>      CMP     R2,R1,#2
>      PGE     R2,TT
>      MOV     R1,#1              // return 1
>      EXIT    R29,R0,#0
>
>      SUB     R30,R1,#2
>      SUB     R1,R1,#1
>      CALL    rfib               // rfib(n-1)
>      MOV     R29,R1
>      MOV     R1,R30
>      CALL    rfib               // rfib(n-2)
>
>      ADD     R1,R1,R29
>      EXIT    R29,R0,#0          // return rfib(n-1)+rfib(n-2)
>

Theoretically, could be done in my ISA something like:
rfib:
ADD -24, SP | MOV LR, R1
MOV 1, R2 | CMPGE 2, R4
BF .L0
MOV.X R12, (SP, 0)
MOV R4, R12 | MOV.Q R1, (SP, 16)
ADD R12, -1, R4
BSR rfib
MOV R2, R13 | ADD R12, -2, R4
BSR rfib
ADD R2, R13, R2
MOV.Q (SP, 16), R1
MOV.X (SP, 0), R12
.L0:
ADD 24, SP
JMP R1

But, my compiler isn't going to pull off anything like this.

Checking, actual BGBCC output:
rfib:
MOV LR, R1
BSR __prolog_0000_0000020000FC
ADD -208, SP
MOV RQ4, RQ13
// tk_shell_chksvar.c:234 { return((x<2)?1:(rfib(x-1)+rfib(x-2))); }
CMPQGE 2, RQ13
BT .L00802A62
MOV 1, RQ10
BRA .L00802A63

..L00802A62:
SUB RQ13, 1, RQ14
MOV RQ14, RQ4
BSR rfib
MOV RQ2, RQ12
SUB RQ13, 2, RQ14
MOV RQ14, RQ4
BSR rfib
MOV RQ2, RQ11
ADD RQ12, RQ11, RQ14
MOV RQ14, RQ10

..L00802A63:
MOV RQ10, RQ2

..L00C108CE:
ADD 208, SP
BRA __epilog_0000_0000020000FC
.balign 4

( Yes, this is more MOV's than "technically necessary", but eliminating
them from the compiler output is "easier said than done" given how some
parts of the compiler work. RQn/RDn are technically equivalent to Rn,
but the compiler distinguishes them as early on the idea was that the
assembler might distinguish instructions based on type, like
EAX/RAX/AX/... on x86-64, but it didn't quite happen this way. )

And, fetching the prolog and epilog:
__prolog_0000_0000020000FC:
ADD -48, SP
MOV.Q R1, (SP, 40)
MOV.X R12, (SP, 16)
MOV.Q R14, (SP, 32)
MOV.X R10, (SP, 0)
RTS

__epilog_0000_0000020000FC:
MOV.Q (SP, 40), R1
MOV.X (SP, 0), R10
MOV.X (SP, 16), R12
MOV.Q (SP, 32), R14
ADD 48, SP
JMP R1

Or, a little closer to the final machine-code (BJX2 Baseline):

rfib: //@03E27C
4911 STC R1, R1
..reloc __prolog_0000_0000020000FC 0F/RELW20_BJX
F000_D000 BSR (PC, 0)
F2F1_D330 ADD -208, SP
18D4 MOV R4, R13
F2DA_C802 CMPQGE 2, R13
..reloc .L00802A62 0F/RELW20_BJX
E000_C000 BT (PC, 0)
DA01 MOV 1, R10
..reloc .L00802A63 0F/RELW20_BJX
F000_C000 BRA (PC, 0)
3000 NOP

..L00802A62: //@03E298
F2ED_11FF ADD R13, -1, R14
F04E_1089 MOV R14, R4
..reloc rfib 0F/RELW20_BJX
F000_D000 BSR (PC, 0)
F4C2_1089 MOV R2, R12 |
F2ED_11FE ADD R13, -2, R14
F04E_1089 MOV R14, R4
..reloc rfib 0F/RELW20_BJX
F000_D000 BSR (PC, 0)
F0B2_1089 MOV R2, R11
F0EC_10B0 ADD R12, R11, R14
F0AE_1089 MOV R14, R10

..L00802A63: //@03E2C0
182A MOV R10, R2

..L00C108CE: //@03E2C2
F2F0_D0D0 ADD 208, SP
..reloc __epilog_0000_0000020000FC 0F/RELW20_BJX
F000_C000 BRA (PC, 0)
3030 BRK

> But then I moved save and restore to the rfib() calls making the last rung
> on the call tree less expensive by 6 memory references each, and by saving
> restoring only 2 registers rather than 3 at the cost of a MOV, so each non-
> leaf level in the call tree saves/restores only 4 registers instead of 6::
>
> rfib:
>      SUB     R2,R1,#2
>      PNLT    R2,TT              // (n-2) < 0
>      MOV     R1,#1              // return 1
>      RET
>
>      ENTER   R30,R0,#0          // just 2 preserved registers.
>      MOV     R30,R2
>      SUB     R1,R1,#1
>      CALL    rfib               // rfib(n-1)
>      MOV     R2,R1
>      MOV     R1,R30
>      MOV     R30,R2
>      CALL    rfib               // rfib(n-2)
>
>      ADD     R1,R1,R30
>      EXIT    R30,R0,#0          // return rfib(n-1)+rfib(n-2)
>
> But then I realized that rfib(n-2) has already been computed by rfib(n-1).
> Since the size of the redundant element on the stack is constant, rfib(n-1)
> makes a location available to rfib(n-2) so only 1 call is required instead
> of 2::
>
> rfib:
>      ENTER   R0,R0,#8
>      STD     #0,[SP+8]        // setup rfib[n-2]
>      CALL    rfib1             // recurse through rfib1
>      EXIT    R0,R0,#8
>
> rfib1:
>      SUB     R2,R1,#2
>      PNLT    R2,TT              // (n-2) < 0
>      MOV     R1,#1              // return 1
>      RET
>
>      ENTER   R0,R0,#8           // no preserved registers just return
> address
>
>      SUB     R1,R1,#1
>      CALL    rfib1              // rfib(n-1)
>      LDD     R2,[SP]            // save rfib(n-2)
>
>      ADD     R1,R1,R30
>      STD     R1,[SP+16]         // restore rfib(n-2)
>      EXIT    R0,R0,#8           // return rfib(n-1)+rfib(n-2)
>
> ....

If you cache the intermediate results of the calculations, this entirely
defeats its use as a benchmark...

Say, because it drops from ~ O(2^n) to ~ O(n)...

Granted, There are plenty of other much more efficient ways of
calculating Fibonacci numbers...

BGB wrote:

> On 12/30/2023 4:05 PM, MitchAlsup wrote:
>> BGB wrote:
>>
>>>    long rfib(long x)
>>>      { return((x<2)?1:(rfib(x-1)+rfib(x-2))); }
>>
>> At first I thought the following looked pretty good:
>>
>> rfib:
>>      ENTER   R29,R0,#0          // just 3 preserved registers.
>>
>>      CMP     R2,R1,#2
>>      PGE     R2,TT
>>      MOV     R1,#1              // return 1
>>      EXIT    R29,R0,#0
>>
>>      SUB     R30,R1,#2
>>      SUB     R1,R1,#1
>>      CALL    rfib               // rfib(n-1)
>>      MOV     R29,R1
>>      MOV     R1,R30
>>      CALL    rfib               // rfib(n-2)
>>
>>      ADD     R1,R1,R29
>>      EXIT    R29,R0,#0          // return rfib(n-1)+rfib(n-2)
>>

> Theoretically, could be done in my ISA something like:
> rfib:
> ADD -24, SP | MOV LR, R1
> MOV 1, R2 | CMPGE 2, R4
> BF .L0
> MOV.X R12, (SP, 0)
> MOV R4, R12 | MOV.Q R1, (SP, 16)
> ADD R12, -1, R4
> BSR rfib
> MOV R2, R13 | ADD R12, -2, R4
> BSR rfib
> ADD R2, R13, R2
> MOV.Q (SP, 16), R1
> MOV.X (SP, 0), R12
> .L0:
> ADD 24, SP
> JMP R1

> But, my compiler isn't going to pull off anything like this.

> Checking, actual BGBCC output:
> rfib:
> MOV LR, R1
> BSR __prolog_0000_0000020000FC
> ADD -208, SP
> MOV RQ4, RQ13
> // tk_shell_chksvar.c:234 { return((x<2)?1:(rfib(x-1)+rfib(x-2))); }
> CMPQGE 2, RQ13
> BT .L00802A62
> MOV 1, RQ10
> BRA .L00802A63

> ..L00802A62:
> SUB RQ13, 1, RQ14
> MOV RQ14, RQ4
> BSR rfib
> MOV RQ2, RQ12
> SUB RQ13, 2, RQ14
> MOV RQ14, RQ4
> BSR rfib
> MOV RQ2, RQ11
> ADD RQ12, RQ11, RQ14
> MOV RQ14, RQ10

> ..L00802A63:
> MOV RQ10, RQ2

> ..L00C108CE:
> ADD 208, SP
> BRA __epilog_0000_0000020000FC
> .balign 4

> ( Yes, this is more MOV's than "technically necessary", but eliminating
> them from the compiler output is "easier said than done" given how some
> parts of the compiler work. RQn/RDn are technically equivalent to Rn,
> but the compiler distinguishes them as early on the idea was that the
> assembler might distinguish instructions based on type, like
> EAX/RAX/AX/... on x86-64, but it didn't quite happen this way. )

> And, fetching the prolog and epilog:
> __prolog_0000_0000020000FC:
> ADD -48, SP
> MOV.Q R1, (SP, 40)
> MOV.X R12, (SP, 16)
> MOV.Q R14, (SP, 32)
> MOV.X R10, (SP, 0)
> RTS

> __epilog_0000_0000020000FC:
> MOV.Q (SP, 40), R1
> MOV.X (SP, 0), R10
> MOV.X (SP, 16), R12
> MOV.Q (SP, 32), R14
> ADD 48, SP
> JMP R1

> Or, a little closer to the final machine-code (BJX2 Baseline):

> rfib: //@03E27C
> 4911 STC R1, R1
> ..reloc __prolog_0000_0000020000FC 0F/RELW20_BJX
> F000_D000 BSR (PC, 0)
> F2F1_D330 ADD -208, SP
> 18D4 MOV R4, R13
> F2DA_C802 CMPQGE 2, R13
> ..reloc .L00802A62 0F/RELW20_BJX
> E000_C000 BT (PC, 0)
> DA01 MOV 1, R10
> ..reloc .L00802A63 0F/RELW20_BJX
> F000_C000 BRA (PC, 0)
> 3000 NOP

> ..L00802A62: //@03E298
> F2ED_11FF ADD R13, -1, R14
> F04E_1089 MOV R14, R4
> ..reloc rfib 0F/RELW20_BJX
> F000_D000 BSR (PC, 0)
> F4C2_1089 MOV R2, R12 |
> F2ED_11FE ADD R13, -2, R14
> F04E_1089 MOV R14, R4
> ..reloc rfib 0F/RELW20_BJX
> F000_D000 BSR (PC, 0)
> F0B2_1089 MOV R2, R11
> F0EC_10B0 ADD R12, R11, R14
> F0AE_1089 MOV R14, R10

> ..L00802A63: //@03E2C0
> 182A MOV R10, R2

> ..L00C108CE: //@03E2C2
> F2F0_D0D0 ADD 208, SP
> ..reloc __epilog_0000_0000020000FC 0F/RELW20_BJX
> F000_C000 BRA (PC, 0)
> 3030 BRK

>> But then I moved save and restore to the rfib() calls making the last rung
>> on the call tree less expensive by 6 memory references each, and by saving
>> restoring only 2 registers rather than 3 at the cost of a MOV, so each non-
>> leaf level in the call tree saves/restores only 4 registers instead of 6::
>>
>> rfib:
>>      SUB     R2,R1,#2
>>      PNLT    R2,TT              // (n-2) < 0
>>      MOV     R1,#1              // return 1
>>      RET
>>
>>      ENTER   R30,R0,#0          // just 2 preserved registers.
>>      MOV     R30,R2
>>      SUB     R1,R1,#1
>>      CALL    rfib               // rfib(n-1)
>>      MOV     R2,R1
>>      MOV     R1,R30
>>      MOV     R30,R2
>>      CALL    rfib               // rfib(n-2)
>>
>>      ADD     R1,R1,R30
>>      EXIT    R30,R0,#0          // return rfib(n-1)+rfib(n-2)
>>
>> But then I realized that rfib(n-2) has already been computed by rfib(n-1).
>> Since the size of the redundant element on the stack is constant, rfib(n-1)
>> makes a location available to rfib(n-2) so only 1 call is required instead
>> of 2::
>>
>> rfib:
>>      ENTER   R0,R0,#8
>>      STD     #0,[SP+8]        // setup rfib[n-2]
>>      CALL    rfib1             // recurse through rfib1
>>      EXIT    R0,R0,#8
>>
>> rfib1:
>>      SUB     R2,R1,#2
>>      PNLT    R2,TT              // (n-2) < 0
>>      MOV     R1,#1              // return 1
>>      RET
>>
>>      ENTER   R0,R0,#8           // no preserved registers just return
>> address
>>
>>      SUB     R1,R1,#1
>>      CALL    rfib1              // rfib(n-1)
>>      LDD     R2,[SP]            // save rfib(n-2)
>>
>>      ADD     R1,R1,R30
>>      STD     R1,[SP+16]         // restore rfib(n-2)
>>      EXIT    R0,R0,#8           // return rfib(n-1)+rfib(n-2)
>>
>> ....

> If you cache the intermediate results of the calculations, this entirely
> defeats its use as a benchmark...

Exactly:: but there is talk of compilers figuring out that rfib() is a
pure function and then the compiler itself performs that transformation.

{{Oh and BTW: Fibonacci should have never been a benchmark (post 1985)}}

> Say, because it drops from ~ O(2^n) to ~ O(n)...

It also reduces the constant term in BigO ! since the code path is shorter.

Back in 1982±, S.E.L. got their compiler to recognize the transcendental
loop as constant, and this one transformation improved Whetstone by 3×.

> Granted, There are plenty of other much more efficient ways of
> calculating Fibonacci numbers...

Make a file as big as you like, and simply index.

BigO( const )

On 30/12/2023 23:19, BGB wrote:

> Granted, There are plenty of other much more efficient ways of
> calculating Fibonacci numbers...
>

Yes. In order of effort.

1) Memoize the recursive function.
2) Return values for (n-1,n-2) as a tuple, and make it single (tail)
recursive. Which a smart compiler could theoretically convert to an
iterative loop.
3) Write an iterative loop yourself.
4) Remember college maths, and solve the recurrence relation to give a
very simple closed form solution.

In fact, given that fib grows exponentially, so we only ever need to
calculate it for small values of n, the only “bad” solution is the naive
double recursive one.

MitchAlsup wrote:
> BGB wrote:
>
>>    long rfib(long x)
>>      { return((x<2)?1:(rfib(x-1)+rfib(x-2))); }
>
> At first I thought the following looked pretty good:
>
> rfib:
>      ENTER   R29,R0,#0          // just 3 preserved registers.
>
>      CMP     R2,R1,#2
>      PGE     R2,TT
>      MOV     R1,#1              // return 1
>      EXIT    R29,R0,#0
>
>      SUB     R30,R1,#2
>      SUB     R1,R1,#1
>      CALL    rfib               // rfib(n-1)
>      MOV     R29,R1
>      MOV     R1,R30
>      CALL    rfib               // rfib(n-2)
>
>      ADD     R1,R1,R29
>      EXIT    R29,R0,#0          // return rfib(n-1)+rfib(n-2)
>
> But then I moved save and restore to the rfib() calls making the last rung
> on the call tree less expensive by 6 memory references each, and by saving
> restoring only 2 registers rather than 3 at the cost of a MOV, so each non-
> leaf level in the call tree saves/restores only 4 registers instead of 6::
>
> rfib:
>      SUB     R2,R1,#2
>      PNLT    R2,TT              // (n-2) < 0
>      MOV     R1,#1              // return 1
>      RET
>
>      ENTER   R30,R0,#0          // just 2 preserved registers.
>      MOV     R30,R2
>      SUB     R1,R1,#1
>      CALL    rfib               // rfib(n-1)
>      MOV     R2,R1
>      MOV     R1,R30
>      MOV     R30,R2
>      CALL    rfib               // rfib(n-2)
>
>      ADD     R1,R1,R30
>      EXIT    R30,R0,#0          // return rfib(n-1)+rfib(n-2)
>
> But then I realized that rfib(n-2) has already been computed by rfib(n-1).
> Since the size of the redundant element on the stack is constant, rfib(n-1)
> makes a location available to rfib(n-2) so only 1 call is required instead
> of 2::

That is the huge/obvious/sneaky (take your pick!) Fib optimization,
since you are effectively getting rid of the O(2^n) exponential
recursion fanout, making it O(n) instead.

The step from that to tail recursion elimination is the only thing
remaining, right?

What you've done is similar to converting

f_n = fib(n)

into

(f_n,f_n_minus_1) = fib2(n)
{ if (n < 2) return (1,1);
(f2,f1) = fib2(n-1);
return (f1+f2,f1);
}

Terje

--
- <Terje.Mathisen at tmsw.no>
"almost all programming can be viewed as an exercise in caching"

On 31/12/2023 09:12, Pancho wrote:
> In fact, given that fib grows exponentially, so we only ever need to
> calculate it for small values of n, the only “bad” solution is the naive
> double recursive one.

Well, I might disagree.

Given that we only ever need to calculate it for small values of n (as
you say, and if that's true) its efficiency doesn't matter. If you need
to compute it often enough that efficiency might become relevant just
put the results in a table and use that. Given that fib goes over 10^12
after about 60 terms it needn't be a very big table.

The naive doubly recursive solution has the huge merit of being
human-readable.

Unlike, say:

double fib( int n )
{
return (pow( (1+sqrt(5.0))/2, n+1 )
- pow( (1-sqrt(5.0))/2, n+1 ))
/ sqrt(5.0);
}

.... which is great for fibs of big values, but not fast for small ones.

--
Cheers,
Daniel.

On 31/12/2023 15:45, Daniel James wrote:
> On 31/12/2023 09:12, Pancho wrote:
>> In fact, given that fib grows exponentially, so we only ever need to
>> calculate it for small values of n, the only “bad” solution is the
>> naive double recursive one.
>
> Well, I might disagree.
>
> Given that we only ever need to calculate it for small values of n (as
> you say, and if that's true) its efficiency doesn't matter. If you need
> to compute it often enough that efficiency might become relevant just
> put the results in a table and use that. Given that fib goes over 10^12
> after about 60 terms it needn't be a very big table.
>

Memoization, is effectively a lookup table, but without the thinking. In
my case, avoiding thinking massively improves code reliability.

> The naive doubly recursive solution has the huge merit of being
> human-readable.
>Yes, but you'd be surprised how quickly double recursion goes bad. I
haven't tested recently, but off the top of my head, the maximum
function call depth would be something like 300, so fib(10) would be a
stack overflow. fib(60) would need a stack function call depth ~ 1e18.

> Unlike, say:
>
>   double fib( int n )
>   {
>       return (pow( (1+sqrt(5.0))/2, n+1 )
>               - pow( (1-sqrt(5.0))/2, n+1 ))
>           / sqrt(5.0);
>   }
>
> ... which is great for fibs of big values, but not fast for small ones.
>

Why isn't it fast for small fibs?

I seem to recall, if statements are more costly than floating point
function calls, like exponential.

On 29/12/2023 22:04, Andy Burns wrote:
> and (perversely?) with -Ofast

Thank you for those. It looks as though there are some other effects
going on apart from the cache size - but I can clearly see the drop when
you run out of cache. Within the caches? Too much noise for me to
understand!

Andy

Pancho wrote:
> On 31/12/2023 15:45, Daniel James wrote:
>> On 31/12/2023 09:12, Pancho wrote:
>>> In fact, given that fib grows exponentially, so we only ever need to
>>> calculate it for small values of n, the only “bad” solution is the
>>> naive double recursive one.
>>
>> Well, I might disagree.
>>
>> Given that we only ever need to calculate it for small values of n (as
>> you say, and if that's true) its efficiency doesn't matter. If you
>> need to compute it often enough that efficiency might become relevant
>> just put the results in a table and use that. Given that fib goes over
>> 10^12 after about 60 terms it needn't be a very big table.
>>
>
> Memoization, is effectively a lookup table, but without the thinking. In
> my case, avoiding thinking massively improves code reliability.
>
>> The naive doubly recursive solution has the huge merit of being
>> human-readable.
>> Yes, but you'd be surprised how quickly double recursion goes bad. I
> haven't tested recently, but off the top of my head, the maximum
> function call depth would be something like 300, so fib(10) would be a
> stack overflow. fib(60) would need a stack function call depth ~ 1e18.
>
>> Unlike, say:
>>
>>    double fib( int n )
>>    {
>>        return (pow( (1+sqrt(5.0))/2, n+1 )
>>                - pow( (1-sqrt(5.0))/2, n+1 ))
>>            / sqrt(5.0);
>>    }
>>
>> ... which is great for fibs of big values, but not fast for small ones.
>>
>
> Why isn't it fast for small fibs?
>
> I seem to recall, if statements are more costly than floating point
> function calls, like exponential.

The recursive fob() function has just a single if () statement, so worst
case you'll get one or two branch mispredicts, while the 2^n recursive
calls (for fib(10 or 11)) will probably each take about as long.

An iterative O(n) version will take 2-4 clcok cycles per iteration, so
for n less than about 40, you'd only have ~100 clock cycles to evaluate
the two pow() calls. (I'm assuming the compiler or programmer have
pre-calculated both (1+sqrt(5))*0.5 and (1-sqrt(5))*0.5 as well as
1/sqrt(5), so that only the two pow() calls remain!)

You need Mitch's ~20-cycle pow() to win this one at one or two-digit N.

Terje

--
- <Terje.Mathisen at tmsw.no>
"almost all programming can be viewed as an exercise in caching"

On 31/12/2023 16:40, Pancho wrote:
> I
> haven't tested recently, but off the top of my head, the maximum
> function call depth would be something like 300, so fib(10) would be a
> stack overflow. fib(60) would need a stack function call depth ~ 1e18.
>

Sorry this was a mistake, the function call depth is only ~n not 2^n.
For fib(n), there are ~2^n function calls, but not 2^n deep on the stack.

On 31/12/2023 17:03, Terje Mathisen wrote:

>
> The recursive fob() function has just a single if () statement, so worst
> case you'll get one or two branch mispredicts, while the 2^n recursive
> calls (for fib(10 or 11)) will probably each take about as long.
>

How many cycles for a branch mispredict? From real life, I have seen if
statements be more expensve than pow, but it was a long time ago, and
possibly code specific. It surprised me at the time, but ever since I
tend to ignore the cost of functions like pow().

> An iterative O(n) version will take 2-4 clcok cycles per iteration, so
> for n less than about 40, you'd only have ~100 clock cycles to evaluate
> the two pow() calls. (I'm assuming the compiler or programmer have
> pre-calculated both (1+sqrt(5))*0.5 and (1-sqrt(5))*0.5 as well as
> 1/sqrt(5), so that only the two pow() calls remain!)
>
> You need Mitch's ~20-cycle pow() to win this one at one or two-digit N.
>

I generally don't worry about optimizing stuff that is very fast.
Premature optimisation is the root of all evil and what not. If you are
really calling it a huge number of times, a lookup is quickest.

On 31/12/2023 16:40, Pancho wrote:
>> The naive doubly recursive solution has the huge merit of being
>> human-readable.
>
> Yes, but you'd be surprised how quickly double recursion goes bad. I
> haven't tested recently, but off the top of my head, the maximum
> function call depth would be something like 300, so fib(10) would be
> a stack overflow. fib(60) would need a stack function call depth ~
> 1e18.

On 64-bit gcc here on Debian (so, 64-bit ints) I can call fib(50)
without blowing the stack, but it takes minutes (on a Ryzen 7).

However, the (signed) 64-bit integer result overflows for fib(n) where
n>45, which is more of a problem.

>> Unlike, say:
>>
>> double fib( int n )
>> {
>> return (pow( (1+sqrt(5.0))/2, n+1 )
>> - pow( (1-sqrt(5.0))/2, n+1 ))
>> / sqrt(5.0);
>> }
>>
>> ... which is great for fibs of big values, but not fast for small
>> ones.
>>
>
> Why isn't it fast for small fibs?

Gut feeling. I haven't timed it, and it may be.

I mean ... obviously it will be as fast for small ones as for large, but
for small ones it won't beat an iterative calculation and for very small
ones I doubt it will beat the naive recursive method.

My point, though, was that clear code is often more important than fast
execution, and while the naive recursive method for fibs gets horribly
slow for large values even it is acceptable for small ones.

> I seem to recall, if statements are more costly than floating point
> function calls, like exponential.

if statements are costly if they invalidate instruction prefetch and/or
jump to code that isn't in the cache ... but optimizers and CPUs are
pretty good at branch prediction and cacheing.

--
Cheers,
Daniel.

Daniel James wrote:

> Unlike, say:

> double fib( int n )
> {
> return (pow( (1+sqrt(5.0))/2, n+1 )
> - pow( (1-sqrt(5.0))/2, n+1 ))
> / sqrt(5.0);
> }

Looks fast to me::

fib:
ADD R2,R1,#1
POWN R3,#(1+sqrt(5.0))/2,R2
POWN R4,#(1-sqrt(5.0))/2,R2
FADD R2,R3,R4
FDIV R1,R2,#sqrt(5)
RET

Terje Mathisen wrote:

> Pancho wrote:
>>> Unlike, say:
>>>
>>>    double fib( int n )
>>>    {
>>>        return (pow( (1+sqrt(5.0))/2, n+1 )
>>>                - pow( (1-sqrt(5.0))/2, n+1 ))
>>>            / sqrt(5.0);
>>>    }
>>>
>>> ... which is great for fibs of big values, but not fast for small ones.
>>>
>>
>> Why isn't it fast for small fibs?
>>
>> I seem to recall, if statements are more costly than floating point
>> function calls, like exponential.

> The recursive fob() function has just a single if () statement, so worst
> case you'll get one or two branch mispredicts, while the 2^n recursive
> calls (for fib(10 or 11)) will probably each take about as long.

> An iterative O(n) version will take 2-4 clcok cycles per iteration, so
> for n less than about 40, you'd only have ~100 clock cycles to evaluate
> the two pow() calls.

POW is 39 cycles (unpipelined) in My 66000 ISA.

> (I'm assuming the compiler or programmer have
> pre-calculated both (1+sqrt(5))*0.5 and (1-sqrt(5))*0.5 as well as
> 1/sqrt(5), so that only the two pow() calls remain!)

Plus a FSUB and FDIV.

> You need Mitch's ~20-cycle pow() to win this one at one or two-digit N.

The Ln2() is 14, the exp2() is another 14, then there is the 90-bit×64-bit
multiply 5-cycles and pipeline latency of 5 cycles (and I seem to be missing
a cycle)

> Terje

Pancho wrote:

> On 31/12/2023 17:03, Terje Mathisen wrote:

>>
>> The recursive fob() function has just a single if () statement, so worst
>> case you'll get one or two branch mispredicts, while the 2^n recursive
>> calls (for fib(10 or 11)) will probably each take about as long.
>>

> How many cycles for a branch mispredict? From real life, I have seen if
> statements be more expensve than pow, but it was a long time ago, and
> possibly code specific. It surprised me at the time, but ever since I
> tend to ignore the cost of functions like pow().

Long latency branch mispredicts are often in the 13-cycle range. So if you
are issuing 2 instructions per cycle into the pipeline, you may be throwing
~40 instructions per mispredict. For example LD-CMP-BR where the LD misses
L1 but hits in L2.

Pancho <Pancho.Jones@proton.me> writes:

> On 30/12/2023 23:19, BGB wrote:
>
>> Granted, There are plenty of other much more efficient ways of
>> calculating Fibonacci numbers...
>
> Yes. In order of effort.
>
> 1) Memoize the recursive function.
> 2) Return values for (n-1,n-2) as a tuple, and make it single (tail)
> recursive. Which a smart compiler could theoretically convert to an
> iterative loop.
> 3) Write an iterative loop yourself.
> 4) Remember college maths, and solve the recurrence relation to give a
> very simple closed form solution.

All of these approaches have drawbacks of one sort or another.

For (1), we end up writing more code and using more memory.
Some problems are a good fit to memoization, but computing
Fibonacci numbers isn't one of the.

For (2), compilers are not so good at dealing with tail
recursion when the return value is a pair of results
rather than a single result.

For (3), iterative code takes more mental effort to see
that it works than does a simple recursive version.

For (4), limitations of floating point cause wrong values to
be produced, for just for 64-bit inputs.

Also, simply storing values in a table takes a fair amount
of memory: almost 3k bytes for 128-bit inputs, or about
3/4 k bytes for 64-bit inputs.

Here is a recursive formulation that is better than the
suggestion given in (2):

typedef unsigned long long NN;

static NN fib2( NN a, NN b, NN k );

NN
fibonacci( NN n ){
return fib2( 1, 0, n );
}

static NN
fib2( NN a, NN b, NN k ){
return k == 0 ? b : fib2( b, a+b, k-1 );
}

At optimization -Os or better, the resulting code looks
very much the same as if a single iterative loop had
been written instead.

Actually we can do much better, but we should do at least
as well as the code shown above.

Daniel James <daniel@me.invalid> writes:

> However, the (signed) 64-bit integer result overflows for fib(n) where
> n>45, which is more of a problem.

fib( 92 ) is still less than 1UL << 63 (or fib(93) for a 64-bit
unsigned type).

And fib( 186 ) still fits in a 128-bit unsigned type.