Fibonacci matrix-exponentiation: Difference between revisions

The '''Fibonacci sequence''' defined with [[wp:Fibonacci_number#Matrix_form|matrix-exponentiation]]:

 

;Task:

 

 

;Extra:

Generate Fibonacci(2¹⁶), Fibonacci(2³²) and Fibonacci(2⁶⁴) using the same method or another one.

 

;Related tasks:

*   [[Fibonacci sequence]]

*   [[Matrix-exponentiation operator]]

 

=={{header|C sharp|C#}}==

===Matrix exponentiation===

using System.IO;

using System.Numerics;

throw new ArgumentException("Illegal matrix dimensions for multiplication.");

}

for (int i = 0; i < A.GetLength(0); ++i)

{

for (int j = 0; j < B.GetLength(1); ++j)

{

{

C[i, j] += A[i, k] * B[k, j];

}

 

{{out}}

<pre>

I have not attempted to calculate the (2^64)th Fibonacci number which appears to be well out of reach using this approach.

{{libheader|GMP(Go wrapper)}}

 

import (

 

fmt.Printf("Took %s\n\n", time.Since(start))

 

{{out}}

<br>

{{trans|Julia}}

 

import (

 

fmt.Printf("Took %s\n\n", time.Since(start))

 

{{out}}

{{trans|Sidef}}

{{trans|Julia}}

 

import (

s1 := firstFibDigits(n, nd2)

if len(s1) < 20 {

} else {

fmt.Printf(" First 20 : %s\n", s1)

}

fmt.Printf("\nTook %s\n", time.Since(start))

 

{{out}}

<pre>

The digits of the 10th Fibonacci number are:

 

The digits of the 100th Fibonacci number are:

=={{header|Haskell}}==

===Matrix exponentiation===

 

seeFib :: Integer -> String

x = fib (2^n)

xs = take 20 $ show x

fib :: Integer -> Integer

fib 0 = 0 -- this line is necessary because "something ^ 0" returns "fromInteger 1", which unfortunately

-- in our case is not our multiplicative identity (the identity matrix) but just a 1x1 matrix of 1

fib n = (last . head . unMat) (Mat [[1, 1], [1, 0]] ^ n)

mult :: Num a => [[a]] -> [[a]] -> [[a]]

mult uss vss = map ((\xs -> if null xs then [] else foldl1 (zipWith (+)) xs) . zipWith (flip (map . (*))) vss) uss

fromInteger n = Mat [[fromInteger n]]

abs = undefined

{{out}}

Timing is on Intel Core i5-4300U CPU, Windows 10 Professional, using GHCi Version 8.6.5:

<pre>

Fibonacci(2^16) = 73199214460290552832 ... 97270190955307463227

Fibonacci(2^32) = 61999319689381859818 ... 39623735538208076347

</pre>

===Matrix exponentiation - printing alternative ===

 

 

seeFib :: Integer -> String

where

p = 2^n

 

fib :: Integer -> Integer

fib 0 = 0 -- this line is necessary because "something ^ 0" returns "fromInteger 1", which unfortunately

-- in our case is not our multiplicative identity (the identity matrix) but just a 1x1 matrix of 1

fib n = (last . head . unMat) (Mat [[1, 1], [1, 0]] ^ n)

mult :: Num a => [[a]] -> [[a]] -> [[a]]

mult uss vss = map ((\xs -> if null xs then [] else foldl1 (zipWith (+)) xs) . zipWith (flip (map . (*))) vss) uss

fromInteger n = Mat [[fromInteger n]]

abs = undefined

{{out}}

Timing is on Intel Core i5-4300U CPU, Windows 10 Professional, using GHCi Version 8.6.5:

<pre>

</pre>

===Matrix exponentiation for a symmetric matrix===

We will use a property of symmetric matrices which commute.

 

At each step of the exponentiation of a symmetric matric, we multiply 2 symmetric matrices which commute.

import System.CPUTime (getCPUTime)

 

 

seeFib :: Integer -> String

where

start = "Fibonacci(2^" ++ (show n) ++") = "

p = 2^n

 

fib :: Integer -> Integer

fib n = zeroOne (power (Mat 1 1 1 0) n)

 

-- for a symmetric matrix

square :: Num a => (Mat a) -> (Mat a)

where

 

-- for 2 symmetric matrices which commute

mult :: Num a => (Mat a) -> (Mat a) -> (Mat a)

 

power :: Num a => (Mat a) -> Integer -> (Mat a)

power _ 0 = Mat 1 0 0 1

power m 1 = m

{{out}}

Timing is on Intel Core i5-4300U CPU, Windows 10 Professional, using GHCi Version 8.6.5:

 

=={{header|Java}}==

Implemented fib and fibMod.

 

import java.math.BigInteger;

import java.util.Arrays;

 

}

 

{{out}}

fib(1,000,000) = 19532821287077577316 ... 68996526838242546875

fib(10,000,000) = 11298343782253997603 ... 86998673686380546875

</pre>

 

the 2^64-th fibonacchi number, due to BigInt overflow. The Binet method actually overflows even with the 2^32-nd fibonacchi number, so the

Lucas method is used as the alternative method.

const b = [big"1" 1; 1 0]

matrixfibonacci(n) = n == 0 ? 0 : n == 1 ? 1 : (b^(n+1))[2,2]

rpad(firstlast(firstbinet(32) * lastfibmod(32)), 45))

println("2^64 ", " "^90, rpad(firstlast(firstbinet(64) * lastfibmod(64)), 45))

<pre>

N Matrix Lucas Mod

2^64 11175807536929528424...17529800348089840187

</pre>

 

=={{header|Perl}}==

{{trans|Sidef}}

use warnings;

 

my $last20 = nth_fib_last_k_digits(2**$n, 20);

printf "F(2^$n) = %s ... %s\n", $first20, $last20;

{{out}}

<pre>F(2^16) = 73199214460290552832 ... 97270190955307463227

F(2^64) = 11175807536929528424 ... 17529800348089840187</pre>

 

Following the general approach of Sidef, and relying on Perl for <code>fibmod</code> function. Borrowed/adapted routines from [http://rosettacode.org/wiki/Ramanujan%27s_constant Ramanujan's_constant] task to allow <code>FatRat</code> calculations throughout. Does not quite meet task spec, as n^64 results not working yet.

use Inline::Perl5;

my $p5 = Inline::Perl5.new();

my $last20 = nth_fib_last_k_digits(2**$n, 20);

printf "F(2^$n) = %s ... %s\n", $first20, $last20

{{out}}

<pre>F(2^16) = 73199214460290552832 ... 97270190955307463227

F(2^16) = 73199214460290552832 ... 97270190955307463227

F(2^32) = 61999319689381859818 ... 39623735538208076347</pre>

 

=={{header|Sidef}}==

Computing the n-th Fibonacci number, using matrix-exponentiation (this function is also built-in):

([[1,1],[1,0]]**n)[0][1]

}

 

 

First and last 20 digits of the n-th Fibonacci number:

 

var f = fibonacci(2**n)

say "F(2^#{n}) = #{a} ... #{b}"

}

{{out}}

<pre>

 

More efficient approach, using Binet's formula for computing the first k digits, combined with the built-in method '''fibmod(n,m)''' for computing the last k digits:

const PHI = (1.25.sqrt + 0.5)

const IHP = -(1.25.sqrt - 0.5)

var last20 = nth_fib_last_k_digits(2**n, 20)

say "F(2^#{n}) = #{first20} ... #{last20}"

{{out}}

<pre>

F(2^32) = 61999319689381859818 ... 39623735538208076347

F(2^64) = 11175807536929528424 ... 17529800348089840187

</pre>