[seqfan] Re: Jacobsthal function

[hv at crypt.org]

Charles Greathouse <charles.greathouse at case.edu> wrote:
:The Jacobsthal function, A048669, is an important sequence which
:appears in many contexts (e.g., bounding the size of prime gaps;
:finding the least prime in an arithmetic progression; bounding the
:state complexity of regular languages). It is defined as the maximal
:distance between numbers relatively prime to n.
:
:How is it computed? The sequence is marked easy, and yet the only
:program listed is very slow, doing calculations on all numbers up to
:n.
:
:A search reveals a body of literature on computing A048670 =
:a(prime(n)#) efficiently, but nothing on the sequence itself. Perhaps
:it is obvious but I can't immediately see it.

I've written some code for this, available on github at
<https://github.com/hvds/seq/tree/master/jacobsthal>.

I also used calculation of prime(n)# to benchmark what I was doing,
but the code itself accepts an arbitrary list of prime factors.

There are two versions of the code there, one in perl (requires a perl
interpreter to run), the other in C (requires a C compiler to compile).
They are both used the same way, accepting an optional '-k<min>' for the
minimum length to start searching at, followed by a list of prime
factors of the n to consider. The C version is about 100 times faster.

Example:
% ./jacobsthal 2 3 5 7 23 31 37
[0.00] Jd = 7 [0.*.*.*.*.*]
[0.00] Jd = 8 [0.*.*.*.*.*]
[0.00] Jd = 9 [0.*.1.*.*.*]
[0.00] Jd = 10 [0.*.1.*.*.*]
[0.00] Jd = 11 [0.0.1.*.*.*]
[0.00] Jd = 23
% ./jacobsthal -k10 2 3 5 7 23 31 37 
[0.00] Jd = 10 [0.*.1.*.*.*]
[0.00] Jd = 11 [0.0.1.*.*.*]
[0.00] Jd = 23
% 

Notes:
1) the result reported is jacobsthal(n)-1, the length of the longest
run of consecutive integers with a common factor to n;
2) '2' is treated separately: the longest run for the odd factors is
considered, then one more than twice the result reported if a factor
of 2 was requested;
3) the number in brackets at the start of each line is the processor
usage at that point;
4) the numbers in brackets at the end of a line represent the offset
within the run of the first multiple of the respective prime factor,
with '*' meaning "don't care, assign it to any unmatched offset".

So the above finds a run of 11, such that the first is divisible by
3 and 5, and the second is divisibly by 7, with three unassigned
offsets that can have the remaining three odd factors assigned in any
order.

The timings suggest the perl code can calculate Jacobsthal(n) for
n = 2 x 3 x 5 x 7 x ... 41 in less than a second, and the C code
(about 100 times faster) for n = 2 x 3 x ... 53. It is quite hard
for me to get accurate timings on this machine though, so best
possible figures may be both faster and more regular than I show.

Generally, the point of using n=k# as the values to test against is
that these represent the hardest combinations of k prime factors,
so the most stringent test of the speed and accuracy of the code.
It is possible that's also what is happening in the existing
literature - there may be nothing there that actually restricts it
to calculating such values.

Hope this helps,

Hugo