Finding Perrin Pseudo Primes_Part 2

   
Recent changes
Table of contents
Links to this page
FRONT PAGE / INDEX

Subscribe!
@ColinTheMathmo
My latest posts can be found here:
Previous blog posts:
Additionally, some earlier writings:

Finding Perrin Pseudo-primes:
Part 2 - 2015/05/17

 
#!/usr/bin/python

n = 2
t = 10

t_prime = 0.0
t_perrin = 0.0

while clock() < int(argv[1]):

    t_prime -= clock()
    is_prime = is_prime_q(n)
    t_prime += clock()
    t_perrin -= clock()
    is_perrin = is_perrin_q(n)
    t_perrin += clock()
    if is_prime != is_perrin:
        print n, is_prime, is_perrin

    if t < clock():
        print 'Tested to %d in %d secs' % (n,t)
        t += 10

    n += 1

print 'Time in prime :', t_prime
print 'Time in Perrin:', t_perrin

# ================
# Edited output:
#
# Tested to 42763 in 100 secs
#
# Time in prime :  0.11
# Time in Perrin: 99.72
#                 =====
#                 99.83
So now we've got the scaffolding of a program to find these Perrin Pseudo-Primes. Here is the main loop of the code, with some simplistic timing added to it. (Note that this code is incomplete and won't run as is).

The output shows that when given 100 seconds to run it gets as far as n=42763, but more importantly, the timing shows that overwhelmingly it spends its time in the routine to test whether or not a number passes the "Perrin Test." So there are a few things we need to do:

  • Make fewer calls to the Perrin Test
  • Make the Perrin Test run faster.

In this post we'll look at a simple technique to accomplish the first of these.

There's a maxim in programming:

You can't make a program
run faster, you can only
make it do less.

One way of making the program do less is to pre-compute things and then use the result multiple times. Another way is to do a quick test to avoid a longer one. We'll do both of those.

Here's the critical observation:

  • If m divides n,
    • and n divides k(n),
      • then m divides k(n).

The second part of the observation is that the Perrin sequence modulo m has to repeat itself after at most $m^3$ steps. There are only $m^3$ possible triples, and the sequence is reversable, so eventually we get the first triple appearing again. Here are the sequences for 2, 3, and 5, up to the point where they start to repeat

  • mod 2 : 1 0 0 1 0 1 1 1 0 0 ...
    • repeats after 7 terms
  • mod 3 : 0 0 2 0 2 2 2 1 1 0 2 1 2 0 0 2 ...
    • repeats after 13 terms
  • mod 5 : 3 0 2 3 2 0 0 2 0 2 2 2 4 4 1 3 0 4 3 4 2 2 1 4 3 0 2 ...
    • repeats after 24 terms 
 
#!/usr/bin/python

def compute_perrin_number(n):

    if n in [2, 3]: return 0
    a0, a1, a2 = 3, 0, 2
    i = n
    while i > 0:
        a0, a1, a2 = a1, a2, (a0+a1) % n
        i -= 1

    return a0

def make_pre_filter(m):

    rtn = [ 3 % m, 0, 2 % m ]
    while True:
        rtn.append( (rtn[-3] + rtn[-2]) % m )
        if rtn[:3] == rtn[-3:]:
            return rtn[:-3]

pre_filters = {}
pre_filter_elements = range(2,104)
for m in pre_filter_elements:
    pre_filters[m] = make_pre_filter(m)

def is_perrin_q(n):

    for m in pre_filter_elements:
        if n % m == 0:
            pre_filter_to_use = pre_filters[m]
            ndx = n % len(pre_filter_to_use)
            if 0 != pre_filter_to_use[ndx]:
                return False

    return compute_perrin_number(n) == 0
So what we can do is this. For each of several primes (and maybe even some composites) we precompute the Perrin sequence modulo that number. Then when we are looking at a number n we look at all the precomputed sequences that divide n and ask if the sequence has a 0 in the right place. If not, then n will not divide k(n). Of course, if it does have a 0 in the right place for every precomputed m that divides n then we will still have to do the computation, but with any luck that will be a lot less.

So here is the result.

  • Without the sieving:
    • Tested to 45761 in 100 secs
    • Perrin computations: 45760

  • With sieving from m=2 to 103:
    • Tested to 109159 in 100 secs
    • Perrin computations: 19004

As you can see, we get more than twice as far, so we're clearly doing a good thing. It's obvious why when we look at the number of Perrin computations we do. Without the sieving every number has the full calculation, but with sieving only 19 in every 109 gets the calculation, or roughly 17.4%, a saving of 82.6%.

But there's still a problem that we can see clearly when we look at the timing:

  • Time in prime : 0.25 secs
  • Time in Perrin : 99.47 secs

Despite saving over 80% we're still spending all our time in the Perrin test routine. In the next installment we'll see how we can speed that up further by using matrices and an adaptation of the Peasant Multiplication Algorithm.

<<<< Prev <<<<
FindingPerrinPseudoPrimes Part1
:
>>>> Next >>>>
Russian Peasant Multiplication

You should follow me on twitter @ColinTheMathmo

Comments

I've decided no longer to include comments directly via the Disqus (or any other) system. Instead, I'd be more than delighted to get emails from people who wish to make comments or engage in discussion. Comments will then be integrated into the page as and when they are appropriate.

If the number of emails/comments gets too large to handle then I might return to a semi-automated system. That's looking increasingly unlikely.

Contents

 

Links on this page

 
Site hosted by Colin and Rachel Wright:
  • Maths, Design, Juggling, Computing,
  • Embroidery, Proof-reading,
  • and other clever stuff.
Suggest a change ( <-- What does this mean?) / Send me email
Front Page / All pages by date / Site overview / Top of page

Universally Browser Friendly     Quotation from
Tim Berners-Lee 		    Valid HTML 3.2!