How can I calculate OEIS A144311 efficiently?

Question

I'm looking for a way to calculate OEIS A144311 efficiently.

In one sense or another, this series considers the number between "relative" twin primes. What do I mean by this?

Well, the number $77$ is relatively prime to {$2,3,5$} even tho it's not actually prime. We're looking for relative twin primes to this set tho: consecutive odds relatively prime to {$2,3,5$}. For instance, $(77,79)$ would be an example of such. Neither number is divisible by {$2,3,5$}.

$78$ then is the number between these. The next smallest number between a similarly relative twin prime pair would be $90$. $(89,91)$ is the next relative (and real) twin prime pair. That gives us a relative twin prime gap of 12 (being $90-78$), and 11 fails between successes. The question is, what is the largest possible twin prime gap relative to just {$2,3,5$}? This so happens to be it. 11 is the max fails for {$2,3,5$}. The series says that 29 is the max fails for {$2,3,5,7$}, and 41 for {$2,3,5,7,11$}, the first 5 primes.

Now that we've well defined these, what would be the best way to calculate them? Also, what is the 17th number in this series? All methods I've used to try to calculate these have resulted in numbers beyond what can be handled well even in something like Mathematica even by the 10th number in the series.

Answer 1

I see basically two approaches, and a hybrid of them that's more efficient than either.

First, you could try out all combinations of residues other than $-1,1$ that the primes could have at the one end of the sequence, and then see how far you get in each case. The number of combinations would be

$$ \prod_{i=3}^n(p_i-2)\;, $$

i.e. $3\cdot5\cdot9\cdot11\cdot15\cdot17\cdot21\cdot27\cdot29\cdot35\cdot39\cdot41\cdot45\cdot51\cdot57\approx4.6\cdot10^{19}$ to reach $p_{17}$, so this is not feasible.

Second, you could start out without fixing residues for the primes and start moving to the right, and whenever you hit a slot that's not covered yet, you pick a prime whose residue hasn't been fixed yet and try the two possible residues to cover the slot. That would require trying out $2^nn!$ different choices, or $2^{17}\cdot17!\approx4.7\cdot10^{19}$ to reach $p_{17}$, funny enough almost the same number, so this isn't feasible, either.

But note that the first approach is costly for large primes and the second approach is costly for many primes, so we can advantageously combine them by adopting the first approach for the small primes and the second approach for the large primes. That is, we try out all combinations of residues for the primes up to, say, $p_7=17$, and then we move right and fix one of the remaining residues whenever we hit an uncovered slot. To reach $p_{17}$, this requires trying out

$$ 3\cdot5\cdot9\cdot11\cdot15\cdot2^{10}\cdot10!\approx8\cdot10^{13} $$

combinations, a whole lot better and within reach of our electronic friends.

Here's java code that applies this approach. I used it to replicate the sequence up to $p_{14}$, which took $1000$ seconds on my MacBook Pro, so you can expect to find the result for $p_{17}$ in about $5760000$ seconds, or about two months. It's probably not a coincidence that this is on the boundary of what's feasible with today's computers; the calculation of the values up to $p_{16}$ (which would take about $2$ days this way) dates from $2009$, and I suspect that they were calculated in precisely this manner.

Answer 2

Taking @joriki's approach, I've now written a version that runs in Mathematica. I'm including now only the latest and just previous revisions to keep this from getting lengthy. The last revision saw a 5x speedup do to the number between twin primes being divisible by 6, and (in the case of my code) the number $2$ and $4$ more than the twins respectively being $\equiv 3 \pmod 6$. Here's the final:

You can run this for free by going to lab.wolframcloud.com/app and creating a free account.

list3[a_, b_, c_] := Union[
  Range[ChineseRemainder[{3, a}, {6, Prime[c]}], b, 6*Prime[c]],
  Range[ChineseRemainder[{3, a + 2}, {6, Prime[c]}], b, 6*Prime[c]]]

tpgaps7[x_, z_] := Module[{plist, tplist, tuples, k, mgap}, 
  plist = Times @@ Table[Prime[y], {y, 1, x}];
  tplist = Complement[Range[3, plist, 6], Flatten[Table[Union[
    Range[2, plist, Prime[n]], 
    Range[4, plist, Prime[n]]], {n, 3, x}]]];
  tuples = Tuples[Table[Table[{a, plist, c}, {a, 0, Prime[c - 1]}], {c, x + 1, z}]];
  k = Length[tuples];
  mgap = AbsoluteTiming[Max[(*Parallel*)Table[Max[Differences@Complement[
    tplist, Flatten[Union @@ Apply[list3, Part[tuples, j], {1}]]]],
    {j, 1, k}]] - 1];
  {mgap}]

Do[Print[tpgaps7[7, i]], {i, 7, 10}]

(*Times when running in Parallel on my Rasp Pi*)
{{0.064885,7,107}}
{{0.097869,8,149}}
{{1.26251,9,203}}
{{35.9214,10,257}}

(*Times for non-parallel online programming lab*)
{{0.003366,7,107}}
{{0.015974,8,149}}
{{0.498724,9,203}}
{{15.4054,10,257}}

Previous DEBUGGED/OPTIMIZED VERSION

list3[a_,b_,c_]:=Union[Range[a,b,Prime[c]],Range[a+2,b,Prime[c]]]

tpgaps6[x_,z_]:=Module[{plist,tplist,tuples,k,mgap},
  plist=Times@@Table[Prime[y],{y,1,x}];
  tplist=Complement[Range[3,plist,6],
    Flatten[Table[Union[
      Range[2,plist,Prime[n]],
      Range[4,plist,Prime[n]]],
    {n,3,x}]]];
  tuples=Tuples[Table[Table[{a,plist,c},{a,0,Prime[c-1]}],{c,x+1,z}]];
  k=Length[tuples];
  mgap=AbsoluteTiming[Max[ParallelTable[
    Max[Differences@Complement[tplist,
      Flatten[Union@@Apply[list3,Part[tuples,j],{1}]]]],{j,1,k}]]-1];
{z,mgap}]

Do[Print[AbsoluteTiming[tpgaps6[7,i]]],{i,7,10}]

(*Times when running in Parallel on my Rasp Pi*)
{0.554687,{7,{0.053122,107}}}
{0.659064,{8,{0.164598,149}}}
{5.49219,{9,{5.00954,203}}}
{173.238,{10,{172.752,257}}}

(*Times for non-parallel online programming lab*)
{0.078032,{7,{0.00339,107}}}
{0.105818,{8,{0.035588,149}}}
{1.45627,{9,{1.38502,203}}}
{48.1555,{10,{48.0839,257}}}

Answer 3

This can be considered as a kind of covering problem, which can be formulated as an integer linear programming problem. Suppose you want to see if $A144311(n) \ge m$.
We ask to cover the set $S = \{1, \ldots, m\}$, where each prime $p$ of the first $n$ primes will cover those members of $S$ congruent to $a_p \pm 1 \mod p$ for some $a_p \in \{0, \ldots, p-1\}$. Given such a covering, if $T$ is an integer $\equiv -a_p \mod p$ for each of these $p$, then $T+1, \ldots, T+m$ is a sequence of $m$ consecutive integers that satisfies the requirement (the existence of $T$ is guaranteed by the Chinese Remainder Theorem).

In the integer linear programming formulation, we have $0-1$ variables $x_{p,a}$ for each of our primes $p$ and $a = 0 \ldots p-1$, with the constraints

$$ \eqalign{\sum_p \sum_{a \equiv i \pm 1 \mod p} x_{p,a} &\ge 1 \ \text{for}\ i=1 \ldots m\cr \sum_{a=0}^{p-1} x_{p,a} &= 1 \ \text{for each}\ p \cr \text{all} & x_{p,a} \in \{0,1\}}$$

For $n = 12$, it took Cplex on my computer about 14 seconds to show that $m=527$ is feasible, and about 22 seconds to show that $m=528$ is infeasible. However, I expect computing times to rise rather rapidly as $n$ increases, although maybe not as rapidly as for exhaustive search methods.

Another possibility is to attack this with a SAT solver.

Feasible solutions might be found using heuristic search methods such as tabu search or simulated annealing, providing lower bounds. Proving upper bounds is likely to be more difficult.