Thanks for taking a look at this issue. 

A fundamental problem here is that there are far too many combinations (in the OP's application) to store in memory all at once. They need to be generated in "blocks" and written to disk. (In my code comments, you'll see many references to the blocks.) My producer/consumer model was actually in three stages, each with a separate thread:

Generate combinations (producer) => Filter (middleman) => Write to disk (consumer).

Handling things in this blockwise fashion and doing the file I/O is why my code is so lengthy.

I've subsequently improved this (but I haven't posted it) so that the Generate stage is subdivided using the Task model. The Filter stage is trivial to parallelize with Threads:-Map. I don't see any point to parallelizing the Write unless one were using multiple disk drives.

In the improved version, I use combinat:-rankcomb and combinat:-unrankcomb to find the starting place of each Task within the canonical ordering of combinations. That way they can be generated iteratively and not all of once.

I am pretty sure that the reason that your parallel version is faster in total CPU time is not better memory management but rather that generating the combinations iteratively via combinat:-nextcomb is substantially more complicated than generating them recursively. Take a look at combinat:-Combinations:-Successor, which is what nextcomb calls.

@dimasje There is unfortunately no command SumTools:-Expand analogous to IntegrationTools:-Expand, although obviously the analogous operation is valid. I only temporarily exchanged Int and sum for the purpose of using Expand. Obviously the exchange is not valid in general.

@dimasje There is unfortunately no command SumTools:-Expand analogous to IntegrationTools:-Expand, although obviously the analogous operation is valid. I only temporarily exchanged Int and sum for the purpose of using Expand. Obviously the exchange is not valid in general.

@dimasje Now would you please re-enter the corrected expression in ordinary keyboard characters so that I can cut and paste it?

@dimasje The two expressions are not equal---the one with sums is the ln of the product.

Would you please reenter your expression using normal keyboard characters? The "Maple Math" format in MaplePrimes is extremely difficult to read.

Is the orbit length known? Your use of the parameter Nit suggests that you know the orbit length when you call the procedure. If you know the length, then you don't need to check whether it reaches the initial point. And if you don't know the length, then there is no need for the parameter.

@candy898 Yes, the two answers are the same. Mine is collected with respect to y, and yours is collected with respect to x. See ?collect .

@Kaare112 That's a good point. I'm glad that you see it that way. See my Answer below.

I must not be understanding your question correctly. Why not just go to the second line, where you set mu, rho, and g, and just change them?

Are the equations in the Word file in proper Maple syntax, with asterisks used for multiplication and "less than or equal" represented as "<="?

I think that it is legitimate. Advanpix is using hardware quadruple precision (128-bit hardware floats). Maple doesn't have that. The quad precision is being compared against Maple's software floats with Digits set to 34.

@brian bovril Okay, here it is with the requested modification. Every time that it gets to a new longest sequence, it prints it and continues searching for the main sequence.

Note that for AD=210, n=11, the only possibilities are sequences starting with 11.

I did find the smallest example of AD = primorial(11), n=11, with the above code in under an hour. I think they were 7-digit  numbers.

IsPrime:= proc(n) option cache; isprime(n) end proc:

APSearch:= proc(AD::posint, n::posint)
     local SP:= 1, longest:= 0, j, k, S;
     do
          SP:= nextprime(SP);
          S:= SP + k*AD $ k = 1..n-1;
          if andmap(IsPrime, [S]) then
               return [SP, S]
          else
               for j from longest to n-1 do
                    if andmap(IsPrime, [S[1..j]]) then
                         longest:= j+1;
                         print([SP,S[1..j]])
                    else
                         break
                    end if
               end do
          end if           
     end do
end proc:   

APSearch(210,11);
                              [2]
                           [13, 223]
                         [13, 223, 433]
                      [13, 223, 433, 643]
                    [13, 223, 433, 643, 853]
                 [13, 223, 433, 643, 853, 1063]
              [47, 257, 467, 677, 887, 1097, 1307]
          [199, 409, 619, 829, 1039, 1249, 1459, 1669]
       [199, 409, 619, 829, 1039, 1249, 1459, 1669, 1879]
    [199, 409, 619, 829, 1039, 1249, 1459, 1669, 1879, 2089]
Warning,  computation interrupted

@brian bovril Okay, here it is with the requested modification. Every time that it gets to a new longest sequence, it prints it and continues searching for the main sequence.

Note that for AD=210, n=11, the only possibilities are sequences starting with 11.

I did find the smallest example of AD = primorial(11), n=11, with the above code in under an hour. I think they were 7-digit  numbers.

IsPrime:= proc(n) option cache; isprime(n) end proc:

APSearch:= proc(AD::posint, n::posint)
     local SP:= 1, longest:= 0, j, k, S;
     do
          SP:= nextprime(SP);
          S:= SP + k*AD $ k = 1..n-1;
          if andmap(IsPrime, [S]) then
               return [SP, S]
          else
               for j from longest to n-1 do
                    if andmap(IsPrime, [S[1..j]]) then
                         longest:= j+1;
                         print([SP,S[1..j]])
                    else
                         break
                    end if
               end do
          end if           
     end do
end proc:   

APSearch(210,11);
                              [2]
                           [13, 223]
                         [13, 223, 433]
                      [13, 223, 433, 643]
                    [13, 223, 433, 643, 853]
                 [13, 223, 433, 643, 853, 1063]
              [47, 257, 467, 677, 887, 1097, 1307]
          [199, 409, 619, 829, 1039, 1249, 1459, 1669]
       [199, 409, 619, 829, 1039, 1249, 1459, 1669, 1879]
    [199, 409, 619, 829, 1039, 1249, 1459, 1669, 1879, 2089]
Warning,  computation interrupted

alias(x= RootOf(_Z^3+_Z+1)):
Normal(1/x^2) mod 7;