www.digitalmars.com         C & C++   DMDScript  

digitalmars.D.learn - matrix and fibonacci

reply "newcomer[bob]" <bob dontknow.com> writes:
The following is a matrix implementation of the fibonacci
algorithm:

int fib(int n)
{
      int M[2][2] = {{1,0},{0,1}}
      for (int i = 1; i < n; i++)
          M = M * {{1,1},{1,0}}
      return M[0][0];
   }

problem is I don't really understand how matrix multiplication
works so i cannot translate it to an equivalent solution in D.
Thank for you assistance in advance.

newcomer[bob]
Mar 08 2012
next sibling parent reply "newcomer[bob]" <bob dontknow.com> writes:
On Friday, 9 March 2012 at 05:50:03 UTC, newcomer[bob] wrote:
 The following is a matrix implementation of the fibonacci
 algorithm:

 int fib(int n)
 {
      int M[2][2] = {{1,0},{0,1}}
      for (int i = 1; i < n; i++)
          M = M * {{1,1},{1,0}}
      return M[0][0];
   }

 problem is I don't really understand how matrix multiplication
 works so i cannot translate it to an equivalent solution in D.
 Thank for you assistance in advance.

 newcomer[bob]

I attempted the following but it does not work: int fib(int n) ulong[2][2] M = [[1,0], [0,1]]; ulong[2][2] C = [[1,1], [1,0]]; foreach (i; 0 .. n) { M[0][0] = M[0][0] * C[0][0] + M[0][0] * C[0][1]; M[0][1] = M[0][1] * C[0][1] + M[0][1] * C[1][1]; M[1][0] = M[0][1] * C[0][0] + M[1][1] * C[0][1]; M[1][1] = M[1][1] * C[0][1] + M[1][1] * C[1][1]; } return M[0][0]; } any ideas what I'm doing wrong? Thanks, Bob
Mar 09 2012
parent sclytrack <sclytrack hotmail.com> writes:
On 03/09/2012 10:51 AM, newcomer[bob] wrote:
 On Friday, 9 March 2012 at 09:22:47 UTC, newcomer[bob] wrote:
 On Friday, 9 March 2012 at 05:50:03 UTC, newcomer[bob] wrote:
 The following is a matrix implementation of the fibonacci
 algorithm:

 int fib(int n)
 {
 int M[2][2] = {{1,0},{0,1}}
 for (int i = 1; i < n; i++)
 M = M * {{1,1},{1,0}}
 return M[0][0];
 }

 problem is I don't really understand how matrix multiplication
 works so i cannot translate it to an equivalent solution in D.
 Thank for you assistance in advance.

 newcomer[bob]

I attempted the following but it does not work: int fib(int n) long[2][2] M = [[1,0], [0,1]]; ulong[2][2] C = [[1,1], [1,0]]; foreach (i; 0 .. n) { M[0][0] = M[0][0] * C[0][0] + M[0][0] * C[0][1]; M[0][1] = M[0][1] * C[0][1] + M[0][1] * C[1][1]; M[1][0] = M[0][1] * C[0][0] + M[1][1] * C[0][1]; M[1][1] = M[1][1] * C[0][1] + M[1][1] * C[1][1]; } return M[0][0]; } any ideas what I'm doing wrong? Thanks, Bob

Turns out that this this is an algorithm for calculating the powers of two up to 2^63. I still cannot find how to modify it to produce the fibonacci sequence though. Any advice would be appreciated. Thanks, Bob

http://en.wikipedia.org/wiki/Matrix_multiplication#Non-technical_example http://en.wikipedia.org/wiki/Fibonacci_number#Matrix_form I don't understand what you are trying to do above.
Mar 09 2012
prev sibling next sibling parent "newcomer[bob]" <bob dontknow.com> writes:
On Friday, 9 March 2012 at 09:22:47 UTC, newcomer[bob] wrote:
 On Friday, 9 March 2012 at 05:50:03 UTC, newcomer[bob] wrote:
 The following is a matrix implementation of the fibonacci
 algorithm:

 int fib(int n)
 {
     int M[2][2] = {{1,0},{0,1}}
     for (int i = 1; i < n; i++)
         M = M * {{1,1},{1,0}}
     return M[0][0];
  }

 problem is I don't really understand how matrix multiplication
 works so i cannot translate it to an equivalent solution in D.
 Thank for you assistance in advance.

 newcomer[bob]

I attempted the following but it does not work: int fib(int n) long[2][2] M = [[1,0], [0,1]]; ulong[2][2] C = [[1,1], [1,0]]; foreach (i; 0 .. n) { M[0][0] = M[0][0] * C[0][0] + M[0][0] * C[0][1]; M[0][1] = M[0][1] * C[0][1] + M[0][1] * C[1][1]; M[1][0] = M[0][1] * C[0][0] + M[1][1] * C[0][1]; M[1][1] = M[1][1] * C[0][1] + M[1][1] * C[1][1]; } return M[0][0]; } any ideas what I'm doing wrong? Thanks, Bob

Turns out that this this is an algorithm for calculating the powers of two up to 2^63. I still cannot find how to modify it to produce the fibonacci sequence though. Any advice would be appreciated. Thanks, Bob
Mar 09 2012
prev sibling next sibling parent Matthias Walter <xammy xammy.homelinux.net> writes:
On 03/09/2012 10:51 AM, newcomer[bob] wrote:
 On Friday, 9 March 2012 at 09:22:47 UTC, newcomer[bob] wrote:
 On Friday, 9 March 2012 at 05:50:03 UTC, newcomer[bob] wrote:
 The following is a matrix implementation of the fibonacci
 algorithm:

 int fib(int n)
 {
     int M[2][2] = {{1,0},{0,1}}
     for (int i = 1; i < n; i++)
         M = M * {{1,1},{1,0}}
     return M[0][0];
  }

 problem is I don't really understand how matrix multiplication
 works so i cannot translate it to an equivalent solution in D.
 Thank for you assistance in advance.

 newcomer[bob]

I attempted the following but it does not work: int fib(int n) long[2][2] M = [[1,0], [0,1]]; ulong[2][2] C = [[1,1], [1,0]]; foreach (i; 0 .. n) { M[0][0] = M[0][0] * C[0][0] + M[0][0] * C[0][1]; M[0][1] = M[0][1] * C[0][1] + M[0][1] * C[1][1]; M[1][0] = M[0][1] * C[0][0] + M[1][1] * C[0][1]; M[1][1] = M[1][1] * C[0][1] + M[1][1] * C[1][1]; } return M[0][0]; } any ideas what I'm doing wrong? Thanks, Bob

Turns out that this this is an algorithm for calculating the powers of two up to 2^63. I still cannot find how to modify it to produce the fibonacci sequence though. Any advice would be appreciated. Thanks, Bob

The idea is not that bad but it contains a bug: Once you modify M[0][1] in the second step of your loop, its changed content is plugged into the third step. You might simply save the contents of M in 4 additional variables and then use these to fill M again (with the result of M*C). In fact it suffices to store three values as all matrices M (and C) are symmetric, that is, M[0][1] == M[1][0]. But I recommend to do this *after* you made the current version work. Matthias
Mar 09 2012
prev sibling next sibling parent reply Timon Gehr <timon.gehr gmx.ch> writes:
On 03/09/2012 06:50 AM, newcomer[bob] wrote:
 The following is a matrix implementation of the fibonacci algorithm:

 int fib(int n) { int M[2][2] = {{1,0},{0,1}} for (int i = 1; i < n;
 i++) M = M * {{1,1},{1,0}} return M[0][0]; }

 problem is I don't really understand how matrix multiplication works
 so i cannot translate it to an equivalent solution in D. Thank for
 you assistance in advance.

 newcomer[bob]

// simple 2x2 matrix type alias int[2][2] Mat; // 2x2 matrix multiplication Mat matmul(Mat a, Mat b){ return [[a[0][0]*b[0][0]+a[0][1]*b[1][0], a[0][0]*b[0][1]+a[0][1]*b[1][1]], [a[1][0]*b[0][0]+a[1][1]*b[1][0], a[0][1]*b[0][1]+a[1][1]*b[1][1]]]; } // implementation of your algorithm int fib(int n)in{assert(n>=0 && n<46);}body{ int M[2][2] = [[1,0],[0,1]]; int F[2][2] = [[1,1],[1,0]]; foreach (i; 0..n) M = matmul(F,M); return M[0][0]; } // faster way of computing matrix power int fi(int n)in{assert(n>=0 && n<46);}body{ Mat M = [[1,0],[0,1]]; Mat F = [[1,1],[1,0]]; for(;n;n>>=1){ if(n&1) M = matmul(F,M); F = matmul(F,F); } return M[0][0]; } // closed form derived from basis transform to eigenbasis int f(int n)in{assert(n>=0 && n<46);}body{ enum sqrt5=sqrt(5.0); return cast(int)((((1+sqrt5)/2)^^(n+1)-((1-sqrt5)/2)^^(n+1))/sqrt5); }
Mar 09 2012
parent Timon Gehr <timon.gehr gmx.ch> writes:
On 03/10/2012 11:00 PM, newcomer[bob] wrote:
 On Friday, 9 March 2012 at 14:07:10 UTC, Timon Gehr wrote:
 On 03/09/2012 06:50 AM, newcomer[bob] wrote:
 The following is a matrix implementation of the fibonacci algorithm:

 int fib(int n) { int M[2][2] = {{1,0},{0,1}} for (int i = 1; i < n;
 i++) M = M * {{1,1},{1,0}} return M[0][0]; }

 problem is I don't really understand how matrix multiplication works
 so i cannot translate it to an equivalent solution in D. Thank for
 you assistance in advance.

 newcomer[bob]


Thanks very much for the assist. All three of these methods though, seem to have a bug.

I just start the sequence from 1, because a quick glance at your implementation showed that it starts from 1. But it seems like actually it would produce the sequence 1,1,1,2,3,...
 fib(), fi() and f() all produce the
 same incorrect result which for lack of better word I'll call an
 "off by one" bug. A call to any of these functions with any
 integer value between 2 and 44 yields the return value of the
 next higher integer, while 45 overflows.

I don't observe any overflowing. But the assertions are a little too tight, 46 would be the first one that does not work.
 For example, 2 => 2, 10
 => 89, and 15 => 987 which are the actual values for 3, 11 and 16.

 Thanks,
 Bob

Probably this is what you want then: // simple 2x2 matrix type alias int[2][2] Mat; // 2x2 matrix multiplication Mat matmul(Mat a, Mat b){ return [[a[0][0]*b[0][0]+a[0][1]*b[1][0], a[0][0]*b[0][1]+a[0][1]*b[1][1]], [a[1][0]*b[0][0]+a[1][1]*b[1][0], a[0][1]*b[0][1]+a[1][1]*b[1][1]]]; } // implementation of your algorithm int fib(int n)in{assert(n>=0 && n<47);}body{ if(!n) return 0; int M[2][2] = [[1,0],[0,1]]; int F[2][2] = [[1,1],[1,0]]; foreach (i; 1..n) M = matmul(F,M); return M[0][0]; } // faster way of computing matrix power int fi(int n)in{assert(n>=0 && n<47);}body{ if(!n) return 0; Mat M = [[1,0],[0,1]]; Mat F = [[1,1],[1,0]]; for(n--;n;n>>=1){ if(n&1) M = matmul(F,M); F = matmul(F,F); } return M[0][0]; } // closed form derived from basis transform to eigenbasis int f(int n)in{assert(n>=0 && n<47);}body{ enum sqrt5=sqrt(5.0); return cast(int)((((1+sqrt5)/2)^^n-((1-sqrt5)/2)^^n)/sqrt5); }
Mar 10 2012
prev sibling parent "newcomer[bob]" <bob dontknow.com> writes:
On Friday, 9 March 2012 at 14:07:10 UTC, Timon Gehr wrote:
 On 03/09/2012 06:50 AM, newcomer[bob] wrote:
 The following is a matrix implementation of the fibonacci 
 algorithm:

 int fib(int n) { int M[2][2] = {{1,0},{0,1}} for (int i = 1; i 
 < n;
 i++) M = M * {{1,1},{1,0}} return M[0][0]; }

 problem is I don't really understand how matrix multiplication 
 works
 so i cannot translate it to an equivalent solution in D. Thank 
 for
 you assistance in advance.

 newcomer[bob]

// simple 2x2 matrix type alias int[2][2] Mat; // 2x2 matrix multiplication Mat matmul(Mat a, Mat b){ return [[a[0][0]*b[0][0]+a[0][1]*b[1][0], a[0][0]*b[0][1]+a[0][1]*b[1][1]], [a[1][0]*b[0][0]+a[1][1]*b[1][0], a[0][1]*b[0][1]+a[1][1]*b[1][1]]]; } // implementation of your algorithm int fib(int n)in{assert(n>=0 && n<46);}body{ int M[2][2] = [[1,0],[0,1]]; int F[2][2] = [[1,1],[1,0]]; foreach (i; 0..n) M = matmul(F,M); return M[0][0]; } // faster way of computing matrix power int fi(int n)in{assert(n>=0 && n<46);}body{ Mat M = [[1,0],[0,1]]; Mat F = [[1,1],[1,0]]; for(;n;n>>=1){ if(n&1) M = matmul(F,M); F = matmul(F,F); } return M[0][0]; } // closed form derived from basis transform to eigenbasis int f(int n)in{assert(n>=0 && n<46);}body{ enum sqrt5=sqrt(5.0); return cast(int)((((1+sqrt5)/2)^^(n+1)-((1-sqrt5)/2)^^(n+1))/sqrt5); }

Thanks very much for the assist. All three of these methods though, seem to have a bug. fib(), fi() and f() all produce the same incorrect result which for lack of better word I'll call an "off by one" bug. A call to any of these functions with any integer value between 2 and 44 yields the return value of the next higher integer, while 45 overflows. For example, 2 => 2, 10 => 89, and 15 => 987 which are the actual values for 3, 11 and 16. Thanks, Bob
Mar 10 2012