4
$\begingroup$

I've been working on this all night, and I have made this go pretty fast, compared to my first iteration of the program, but now I'm out of ideas.

I'm trying to write a program to test (by good old fashioned brute force) if there is a composite positive integer $n \equiv 3 \pmod {10}$ or $n \equiv 7 \pmod {10}$ where $n$ is a Fermat Pseudoprime, and $n$ divides the $n+1^{\rm th}$ Fibonacci number.

Here is what I have so far:

$HistoryLength = 5; startingIter = 0; iters = startingIter; listLength = 10000; maxIters = 20; nextBlock[start_, listlength_] := Block[ {a, i}, a = DeleteCases[ Table[(3 + 10 i), {i, start, start + listlength - 1}], _?PrimeQ]; a = Join[a, DeleteCases[ Table[(7 + 10 i), {i, start, start + listlength - 1}], _?PrimeQ]]; a = DeleteCases[a, Except[_?(Divisible[Fibonacci[# + 1], #] &)]]; a = DeleteCases[a, Except[_?(Divisible[2^(# - 1) - 1, #] &)]]; Print["Found ", Length[a], " candidates (3 + 10n) or (7 + 10n) between ", TraditionalForm[start], " and ", (start + listLength - 1), "."]; a ]; Monitor[ Do[ Print[ "Iteration: ", ToString[iters + 1], ". Currently testing around: ", ToString[7 + 10*(1 + iters*listLength)], " and a Fibbonaci number with : ", Length[IntegerDigits[Fibonacci[7 + 10*(1 + iters*listLength)]]], " digits: ", nextBlock[1 + iters*listLength, listLength] ], {iters, startingIter, maxIters}], ProgressIndicator[iters, {startingIter, maxIters}] ]

The strategy here is that each time nextBlock[] is called, it first creates a list of all composite integers in the desired range which are congruent to $3$ or $7 \mod 10$. Next, it removes any which do not divide the appropriate Fibonacci number, and then removes any that are not Fermat Pseudo-prime. It's fast enough on my system to do 20 iterations in a reasonable amount of time.

Basically I have no clue what to do next to speed this up. I've looked into Compile[] but the Fibonacci numbers in question have hundreds of thousands of digits, so I can't do that. If anyone is good with parallelization, I could use some guidance since the documentation is difficult to follow. I have an i7 CPU and access to another PC with Mathematica and 4 cores, so if I could figure that out I would probably speed things up a bunch.

I haven't done anything with GPU acceleration in many years, is it possible to port GMP onto CUDA or OpenCL?

Improve this question
$\endgroup$
5
  • $\begingroup$ Since it's a brute force search anyway, you might as well use FindInstance. Could you verify that this is equivalent: FindInstance[ (Mod[n, 10] == 3 || Mod[n, 10] == 7) && Mod[Fibonacci[n + 1], n] == 0 && Mod[2^(n - 1) - 1, n] == 0 && Not[Element[n, Primes]] && n > 0, n, Integers, 1] It will probably sit there for billions of years or forever. Are you sure there is a single number that satisfies these constraints? If you want a faster solution you'll have to use a more intelligent algorithm than brute force. $\endgroup$ Commented May 25, 2015 at 15:33
  • $\begingroup$ It's an open problem. Such a number would be a counterexample to a certain probabilistic Primality test. Really, this is more an interesting exercise for me to get better with Mathematica. $\endgroup$ Commented May 25, 2015 at 15:41
  • $\begingroup$ Well, I just ran FindInstance for a good length of time and it failed to either prove no such number existed or find any number... so my approach doesn't work :( $\endgroup$ Commented May 25, 2015 at 16:05
  • $\begingroup$ It is known that any such number is greater than $2^{64}$ $\endgroup$ Commented May 25, 2015 at 16:08
  • $\begingroup$ ahh, well that explains it : Fibonacci[2^64 + 1] gives Overflow[ ] I'm thinking you might be able to leverage the Fibonacci square test somehow to avoid calculating the Fib altogether, if either (5 (k+1)^2 + 4) or (5 (k+1)^2 - 4) is square, k+1 is a Fibonacci number... tricky $\endgroup$ Commented May 25, 2015 at 16:24

1 Answer 1

Reset to default
8
$\begingroup$
  1. Never use pattern matching unless you absolutely have to. Using Cases instead of Select can make a huge difference.
  2. Vectorize operations. Use Range instead of Table if you can.
  3. Test several things at once. And[test1, test2, test3] will abort when it can for maximum efficiency ("short circuit evaluation").

Taking this into consideration your code looks like this:

nextBlock2[start_, listlength_] := Module[{a}, a = 10 Range[start, start + listlength - 1]; a = Join[3 + a, 7 + a]; Select[a, And[ Not@*PrimeQ, Divisible[Fibonacci[# + 1], #] &, Divisible[2^(# - 1) - 1, #] & ]] ]

And it should be much faster.

Improve this answer
$\endgroup$
2
  • $\begingroup$ This is great, thanks for your help! The only other change I made was to switch the order of Fermat-Pseudoprime testing and fibonacci testing since the fibonacci numbers are so big, so I'm eliminating more with the first test. $\endgroup$ Commented May 25, 2015 at 15:51
  • 1
    $\begingroup$ One more thing: I'd use PowerMod[] instead of using Divisible[] on a big pseudoprime for testing Fermat's little theorem. $\endgroup$ Commented May 25, 2015 at 16:10

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.