Thread overview
matrix and fibonacci
Mar 09, 2012
newcomer[bob]
Mar 09, 2012
newcomer[bob]
Mar 09, 2012
newcomer[bob]
Mar 09, 2012
sclytrack
Mar 09, 2012
Matthias Walter
Mar 09, 2012
Timon Gehr
Mar 10, 2012
newcomer[bob]
Mar 10, 2012
Timon Gehr
March 09, 2012
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]
March 09, 2012
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
March 09, 2012
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
March 09, 2012
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.

March 09, 2012
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
March 09, 2012
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);
}
March 10, 2012
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

March 10, 2012
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);
}