parforstep() available, but parsumstep() is missing

[hermann@stamm-wilbrandt.de]

I would need parsumstep() to improved already impressive parsum() 
runtimes.


I found parforstep() and simulated not yet available functionality. But 
result is poor, although only 20% isprime() tests are done, parforstep 
overhead is unacceptable:

hermann@7950x:~$ numactl -C 0-15 gp -q
? #
   timer = 1 (on)
? parsum(k=1,10^9,isprime(k^2+(k+1)^2))
cpu time = 6min, 49,486 ms, real time = 25,606 ms.
68588950
? c=0;parforstep(k=5,10^9,5,isprime(k^2+(k+1)^2),r,c+=r);c
cpu time = 5min, 40,308 ms, real time = 2min, 39,913 ms.
22858489
?

Is there a better use of parforstep() than what I did?

If not, is implementing and adding parsumstep() an option?


Why do I need that?

parsum(k=1,10^9,isprime(k^2+(k+1)^2)) is fast, but tests too much values 
of k. For values k=1 (mod 5) and k=3 (mod 5) the term k^2+(k+1)^2 is 
divisible by 5, so no isprime() test needed. The equivalent computation 
(+1 needed since 5=1^2+2^2 needs to be counted as prime) would be this, 
and should reduce runtime by 40%:

{
  1
  + parsumstep(k=5,10^9,5,isprime(k^2+(k+1)^2))
  + parsumstep(k=2,10^9,5,isprime(k^2+(k+1)^2))
  + parsumstep(k=4,10^9,5,isprime(k^2+(k+1)^2))
}


Always 2 values do not need to be tested, for 5 as discusssed, and as 
well for 13, 17, 29, 37, 41, 53 listed here with their two values:
https://oeis.org/A012132


With only using the first five primes =1 (mod 4) would have the 
potential to reduce isprime() tests not by 40% but by 60%, although I 
don't know how the parsumstep()s could be combined for that.

? c=0;forprimestep(p=5,40,4,c+=1);c
5
? s=0;forprimestep(p=5,40,4,s+=(2/p-2*s/p));s*1.0
0.60547456490661358815517030527172515002
?


Regards,

Hermann.