Faster way to generate a sequence of numbers

Question

There is a sequence of Numbers that satisfies

(1) $a_1+p\geq 0, a_2+p=0, p\in R$
(2) $a_{4n}>a_{4 n-1}, (n=1,2,...)$
(3) $a_{m+n}\in\{a_m+a_n+p, a_m+a_n+p+1\}, (m,n=1,2...)$

I tried to list the first k items using the code below,

k = 6;
eqn = Table[And@@Function[{m,n},And[a[1]+p>=0,a[2]==-p,a[4]>a[3],a[8]>a[7],
 a[m+n]==a[m]+a[n]+p||a[m+n]==a[m]+a[n]+p+1]]@@@IntegerPartitions[i,{2}],{i,2,k}]
Solve[eqn,Array[a,k],Reals]//Normal

{{a[1]->-p,a[2]->-p,a[3]->-p,a[4]->1-p,a[5]->1-p,a[6]->1-p}}

When I set k to be greater than 6, it has been running for a long time. Do you have a faster way to generate more terms of this sequence?

Answer 1

Here is one suggestion that might work. It seems (on a cursory glance of the conditions making up the sequence) that previous solutions are nested in current ones. For that reason, it probably makes sense to define this recursively. Below, I have worked out a klugey solution based on your current code. It could be made much nicer, but I believe that it works. Part of the kluge is in the following, where in order to make the solutions work, I need to get rid of the conditions a[4] > a[3] and a[8] > a[7] when they are unnecessary:

generateConditions[k_ /; 2 <= k <= 3] := 
 Table[And @@ 
   Function[{m, n}, 
     And[a[1] + p >= 0, a[2] == -p, 
      a[m + n] == a[m] + a[n] + p || 
       a[m + n] == a[m] + a[n] + p + 1]] @@@ 
    IntegerPartitions[i, {2}], {i, 2, k}]
generateConditions[k_ /; 4 <= k <= 7] := 
 Table[And @@ 
   Function[{m, n}, 
     And[a[1] + p >= 0, a[2] == -p, a[4] > a[3], 
      a[m + n] == a[m] + a[n] + p || 
       a[m + n] == a[m] + a[n] + p + 1]] @@@ 
    IntegerPartitions[i, {2}], {i, 2, k}]
generateConditions[k_ /; 8 <= k] := 
 Table[And @@ 
   Function[{m, n}, 
     And[a[1] + p >= 0, a[2] == -p, a[4] > a[3], a[8] > a[7], 
      a[m + n] == a[m] + a[n] + p || 
       a[m + n] == a[m] + a[n] + p + 1]] @@@ 
    IntegerPartitions[i, {2}], {i, 2, k}]

Then, the solutions are recursively defined with the k = 2 case as the base case:

ClearAll@generateSolutions
generateSolutions[2] = Solve[generateConditions[2], Array[a, 2], Reals]

generateSolutions[k_ /; k > 2] := generateSolutions[k] = Module[
   {sols, solOld = generateSolutions[k - 1]}
   , sols = Table[
     Join[sol, #] & /@ 
      Solve[generateConditions[k] /. sol // DeleteDuplicates // 
        Simplify, a[k], Reals]
     , {sol, solOld}]
   ; Cases[sols, {__, a[k] -> _}, Infinity]
   ]

You can run the code to get the following, very quickly. (The k = 7 case is the second one to generate two solutions, and it takes a very long time to solve with the original case, which is why I started there.)

generateSolutions[7]
(* {{a[1] -> -p, a[2] -> -p, a[3] -> -p, a[4] -> 1 - p, a[5] -> 1 - p, a[6] -> 1 - p, a[7] -> 1 - p},
    {a[1] -> -p, a[2] -> -p, a[3] -> -p, a[4] -> 1 - p, a[5] -> 1 - p, a[6] -> 1 - p, a[7] -> 2 - p}} *)

and

generateSolutions[8]
(* {{a[1] -> -p, a[2] -> -p, a[3] -> -p, a[4] -> 1 - p, a[5] -> 1 - p, a[6] -> 1 - p, a[7] -> 1 - p, a[8] -> 2 - p}} *)