Pisano period

From Rosetta Code
Revision as of 11:45, 3 March 2020 by PureFox (talk | contribs) (Added Go)
Pisano period is a draft programming task. It is not yet considered ready to be promoted as a complete task, for reasons that should be found in its talk page.

The Fibonacci sequence taken modulo 2 is a periodic sequence of period 3 : 0, 1, 1, 0, 1, 1, ...

For any integer n, the Fibonacci sequence taken modulo n is periodic and the period is called the Pisano period.

Let call pisano, the Pisano period (pisano(2) = 3). If m and n are coprime, pisano(m*n) = lcm(pisano(m),pisano(n)).

Therefore calculate the Pisano period of an integer m is accomplished by calculating the Pisano periods of the prime powers in the prime decomposition of m.

Task

Write 2 functions: pisanoPrime(p,k) and pisano(m).

pisanoPrime(p,k) return the Pisano period of pk where p is prime and k is a positive integer.

pisano(m) uses pisanoPrime to return the Pisano period of m where m is a positive integer.

Print pisanoPrime(p,2) for every prime lower than 15.

Print pisanoPrime(p,1) for every prime lower than 180.

Print pisano(m) for every integer from 2 to 180.

Related tasks


Go

<lang go>package main

import "fmt"

func gcd(a, b uint) uint {

   if b == 0 {
       return a
   }
   return gcd(b, a%b)

}

func lcm(a, b uint) uint {

   return a / gcd(a, b) * b

}

func ipow(x, p uint) uint {

   prod := uint(1)
   for p > 0 {
       if p&1 != 0 {
           prod *= x
       }
       p >>= 1
       x *= x
   }
   return prod

}

// Gets the prime decomposition of n. func getPrimes(n uint) []uint {

   var primes []uint
   for i := uint(2); i <= n; i++ {
       div := n / i
       mod := n % i
       for mod == 0 {
           primes = append(primes, i)
           n = div
           div = n / i
           mod = n % i
       }
   }
   return primes

}

// Not particularly efficient but suffices here. func isPrime(n uint) bool {

   if n%2 == 0 {
       return n == 2
   }
   if n%3 == 0 {
       return n == 3
   }
   primes := getPrimes(n)
   if len(primes) == 1 && primes[0] == n {
       return true
   }
   return false

}

// Calculates the Pisano period of 'm' from first principles. func pisanoPeriod(m uint) uint {

   var p, c uint = 0, 1
   for i := uint(0); i < m*m; i++ {
       p, c = c, (p+c)%m
       if p == 0 && c == 1 {
           return i + 1
       }
   }
   return 1

}

// Calculates the Pisano period of p^k where 'p' is prime and 'k' is a positive intege. func pisanoPrime(p uint, k uint) uint {

   if !isPrime(p) || k == 0 {
       return 0 // can't do this one
   }
   return pisanoPeriod(ipow(p, k))

}

// Calculates the Pisano period of 'm' using pisanoPrime. func pisano(m uint) uint {

   primes := getPrimes(m)
   primePowers := make(map[uint]uint)
   for _, p := range primes {
       primePowers[p]++
   }
   var pps []uint
   for k, v := range primePowers {
       pps = append(pps, pisanoPrime(k, v))
   }
   if len(pps) == 1 {
       return pps[0]
   }
   f := pps[0]
   for i := 1; i < len(pps); i++ {
       f = lcm(f, pps[i])
   }
   return f

}

func main() {

   for p := uint(2); p < 15; p++ {
       pp := pisanoPrime(p, 2)
       if pp > 0 {
           fmt.Printf("pisanoPrime(%2d: 2) = %d\n", p, pp)
       }
   }
   fmt.Println()
   for p := uint(2); p < 180; p++ {
       pp := pisanoPrime(p, 1)
       if pp > 0 {
           fmt.Printf("pisanoPrime(%3d: 1) = %d\n", p, pp)
       }
   }
   fmt.Println()
   fmt.Println("pisano(n) for integers 'n' from 2 to 180 are:")
   for n := uint(2); n <= 180; n++ {
       fmt.Printf("%3d ", pisano(n))
       if n != 2 && (n-1)%15 == 0 {
           fmt.Println()
       }
   }
   fmt.Println()

}</lang>

Output:
pisanoPrime( 2: 2) = 6
pisanoPrime( 3: 2) = 24
pisanoPrime( 5: 2) = 100
pisanoPrime( 7: 2) = 112
pisanoPrime(11: 2) = 110
pisanoPrime(13: 2) = 364

pisanoPrime(  2: 1) = 3
pisanoPrime(  3: 1) = 8
pisanoPrime(  5: 1) = 20
pisanoPrime(  7: 1) = 16
pisanoPrime( 11: 1) = 10
pisanoPrime( 13: 1) = 28
pisanoPrime( 17: 1) = 36
pisanoPrime( 19: 1) = 18
pisanoPrime( 23: 1) = 48
pisanoPrime( 29: 1) = 14
pisanoPrime( 31: 1) = 30
pisanoPrime( 37: 1) = 76
pisanoPrime( 41: 1) = 40
pisanoPrime( 43: 1) = 88
pisanoPrime( 47: 1) = 32
pisanoPrime( 53: 1) = 108
pisanoPrime( 59: 1) = 58
pisanoPrime( 61: 1) = 60
pisanoPrime( 67: 1) = 136
pisanoPrime( 71: 1) = 70
pisanoPrime( 73: 1) = 148
pisanoPrime( 79: 1) = 78
pisanoPrime( 83: 1) = 168
pisanoPrime( 89: 1) = 44
pisanoPrime( 97: 1) = 196
pisanoPrime(101: 1) = 50
pisanoPrime(103: 1) = 208
pisanoPrime(107: 1) = 72
pisanoPrime(109: 1) = 108
pisanoPrime(113: 1) = 76
pisanoPrime(127: 1) = 256
pisanoPrime(131: 1) = 130
pisanoPrime(137: 1) = 276
pisanoPrime(139: 1) = 46
pisanoPrime(149: 1) = 148
pisanoPrime(151: 1) = 50
pisanoPrime(157: 1) = 316
pisanoPrime(163: 1) = 328
pisanoPrime(167: 1) = 336
pisanoPrime(173: 1) = 348
pisanoPrime(179: 1) = 178

pisano(n) for integers 'n' from 2 to 180 are:
  3   8   6  20  24  16  12  24  60  10  24  28  48  40  24 
 36  24  18  60  16  30  48  24 100  84  72  48  14 120  30 
 48  40  36  80  24  76  18  56  60  40  48  88  30 120  48 
 32  24 112 300  72  84 108  72  20  48  72  42  58 120  60 
 30  48  96 140 120 136  36  48 240  70  24 148 228 200  18 
 80 168  78 120 216 120 168  48 180 264  56  60  44 120 112 
 48 120  96 180  48 196 336 120 300  50  72 208  84  80 108 
 72  72 108  60 152  48  76  72 240  42 168 174 144 120 110 
 60  40  30 500  48 256 192  88 420 130 120 144 408 360  36 
276  48  46 240  32 210 140  24 140 444 112 228 148 600  50 
 36  72 240  60 168 316  78 216 240  48 216 328 120  40 168 
336  48 364 180  72 264 348 168 400 120 232 132 178 120 

Julia

The factorization method gives about a 3X speedup over just looking for the 0, 1 that marks the repeating sequence. <lang julia>using Primes

function pisano(p)

   p < 2 && return 1
   lastn, n = 0, 1
   for i in 0:p^2
       lastn, n = n, (lastn + n) % p
       if lastn == 0 && n == 1
           return i + 1
       end
   end
   return 1

end

pisanoprime(p, k) = (@assert(isprime(p)); pisano(p^k)) primedecomp(n) = [p for p in collect(factor(n))] pisanotask(n) = mapreduce(p -> pisanoprime(p[1], p[2]), lcm, primedecomp(n), init=1)

for i in 1:15

   if isprime(i)
       println("pisanoPrime($i, 2) = ", pisanoprime(i, 2))
   end

end

for i in 1:180

   if isprime(i)
       println("pisanoPrime($i, 1) = ", pisanoprime(i, 1))
   end

end


println("\nPisano(n) for n from 2 to 180:\n", [pisano(i) for i in 2:180]) @time [pisano(i) for i in 9999500:10000000] println("\nPisano(n) using pisanoPrime for n from 2 to 180:\n", [pisanotask(i) for i in 2:180])

@time [pisanotask(i) for i in 9999500:10000000]</lang>

Output:
pisanoPrime(2, 2) = 6
pisanoPrime(3, 2) = 24
pisanoPrime(5, 2) = 100
pisanoPrime(7, 2) = 112
pisanoPrime(11, 2) = 110
pisanoPrime(13, 2) = 364
pisanoPrime(2, 1) = 3
pisanoPrime(3, 1) = 8
pisanoPrime(5, 1) = 20
pisanoPrime(7, 1) = 16
pisanoPrime(11, 1) = 10
pisanoPrime(13, 1) = 28
pisanoPrime(17, 1) = 36
pisanoPrime(19, 1) = 18
pisanoPrime(23, 1) = 48
pisanoPrime(29, 1) = 14
pisanoPrime(31, 1) = 30
pisanoPrime(37, 1) = 76
pisanoPrime(41, 1) = 40
pisanoPrime(43, 1) = 88
pisanoPrime(47, 1) = 32
pisanoPrime(53, 1) = 108
pisanoPrime(59, 1) = 58
pisanoPrime(61, 1) = 60
pisanoPrime(67, 1) = 136
pisanoPrime(71, 1) = 70
pisanoPrime(73, 1) = 148
pisanoPrime(79, 1) = 78
pisanoPrime(83, 1) = 168
pisanoPrime(89, 1) = 44
pisanoPrime(97, 1) = 196
pisanoPrime(101, 1) = 50
pisanoPrime(103, 1) = 208
pisanoPrime(107, 1) = 72
pisanoPrime(109, 1) = 108
pisanoPrime(113, 1) = 76
pisanoPrime(127, 1) = 256
pisanoPrime(131, 1) = 130
pisanoPrime(137, 1) = 276
pisanoPrime(139, 1) = 46
pisanoPrime(149, 1) = 148
pisanoPrime(151, 1) = 50
pisanoPrime(157, 1) = 316
pisanoPrime(163, 1) = 328
pisanoPrime(167, 1) = 336
pisanoPrime(173, 1) = 348
pisanoPrime(179, 1) = 178

Pisano(n) for n from 2 to 180:
[3, 8, 6, 20, 24, 16, 12, 24, 60, 10, 24, 28, 48, 40, 24, 36, 24, 18, 60, 16, 30, 48, 24, 100, 84, 72, 48, 14, 120, 30, 48, 40, 36, 80, 24, 76, 18, 56, 60, 40, 48, 88, 30, 120, 48, 32, 24, 112, 300, 72, 84, 108, 72, 20, 48, 72, 42, 58, 120, 60, 30, 48, 96, 140, 120, 136, 36, 48, 240, 70, 24, 148, 228, 200, 18, 80, 168, 78, 120, 216, 120, 168, 48, 180, 264, 56, 60, 44, 120, 112, 48, 120, 96, 180, 48, 196, 336, 120, 300, 50, 72, 208, 84, 80, 108, 72, 72, 108, 60, 152, 48, 76, 72, 240, 42, 168, 174, 144, 120, 110, 60, 40, 30, 500, 48, 256, 192, 88, 420, 130, 120, 144, 408, 360, 36, 276, 48, 46, 240, 32, 210, 140, 24, 140, 444, 112, 228, 148, 600, 50, 36, 72, 240, 60, 168, 316, 78, 216, 240, 48, 216, 328, 120, 40, 168, 336, 48, 364, 180, 72, 264, 348, 168, 400, 120, 232, 132, 178, 120]
 13.821726 seconds (7 allocations: 4.281 KiB)

Pisano(n) using pisanoPrime for n from 2 to 180:
[3, 8, 6, 20, 24, 16, 12, 24, 60, 10, 24, 28, 48, 40, 24, 36, 24, 18, 60, 16, 30, 48, 24, 100, 84, 72, 48, 14, 120, 30, 48, 40, 36, 80, 24, 76, 18, 56, 60, 40, 48, 88, 30, 120, 48, 32, 24, 112, 300, 72, 84, 108, 72, 20, 48, 72, 42, 58, 120, 60, 30, 48, 96, 140, 120, 136, 36, 48, 240, 70, 24, 148, 228, 200, 18, 80, 168, 78, 120, 216, 120, 168, 48, 180, 264, 56, 60, 44, 120, 112, 48, 120, 96, 180, 48, 196, 336, 120, 300, 50, 72, 208, 84, 80, 108, 72, 72, 108, 60, 152, 48, 76, 72, 240, 42, 168, 174, 144, 120, 110, 60, 40, 30, 500, 48, 256, 192, 88, 420, 130, 120, 144, 408, 360, 36, 276, 48, 46, 240, 32, 210, 140, 24, 140, 444, 112, 228, 148, 600, 50, 36, 72, 240, 60, 168, 316, 78, 216, 240, 48, 216, 328, 120, 40, 168, 336, 48, 364, 180, 72, 264, 348, 168, 400, 120, 232, 132, 178, 120]
  4.208398 seconds (4.23 k allocations: 249.266 KiB)