11149 - Power of Matrix

All about problems in Volume 111. If there is a thread about your problem, please use it. If not, create one with its number in the subject.

Moderator: Board moderators

Darko
Guru
Posts: 580
Joined: Fri Nov 11, 2005 9:34 am
Location: Calgary, Canada

Post by Darko » Tue Jan 02, 2007 8:21 am

This question will probably be silly, but here it goes...

I've never done the matrix exponentiation, so I didn't even try this problem. But what I did was try if what works for integers worked for that matrix in the example. Well, it did work, but I am not sure I could really generalize (because matrix multiplication is not commutative), although it doesn't seem that crazy.

In that example, you can get the answer using the following equation:
x^(k+1) - 1 = (x - 1) * (x^k + ... + x + 1)

Now, it works for 1x1 matrices (heh) and I just tried it for kicks, but it did work for that matrix ({{0,2,0},{0,0,2},{0,0,0}}). To be honest, I was surprised I got the right result. And I am not sure how to deal with a more complicated inverse, this one was nice.

So - can the equation above be generalized for n>1?

Observer
Guru
Posts: 570
Joined: Sat May 10, 2003 4:20 am
Location: Hong Kong

Post by Observer » Tue Jan 02, 2007 8:48 am

Darko wrote:So - can the equation above be generalized for n>1?
The answer is YES.

Recall that for the formula

Code: Select all

                    x^(k+1) - 1
x^k + ... + x + 1 = -----------
                       x - 1
It does not work for x = 1 (in which case denominator = 0).

Now back to the matrix case. The formula below works (you can check it by direct computation):

Code: Select all

(A - I)(A^k + ... + I) = A^(k+1) - I
 
Therefore if A-I is non-singular (i.e. det(A-I) != 0 <=> inverse exists), then

Code: Select all

(A^k + ... + I) = (A - I)^(-1) * (A^(k+1) - I),
                              where (A - I)^(-1) means the inverse of A-I.
But if A-I is singular, then you'll have to think of other methods~ :)
Last edited by Observer on Tue Jan 02, 2007 9:18 am, edited 2 times in total.
7th Contest of Newbies
Date: December 31st, 2011 (Saturday)
Time: 12:00 - 16:00 (UTC)
URL: http://uva.onlinejudge.org

User avatar
cytmike
Learning poster
Posts: 95
Joined: Mon Apr 26, 2004 1:23 pm
Location: Hong Kong and United States
Contact:

Post by cytmike » Tue Jan 02, 2007 8:57 am

removed
Last edited by cytmike on Wed Jan 17, 2007 1:53 am, edited 1 time in total.
Impossible is Nothing.

Darko
Guru
Posts: 580
Joined: Fri Nov 11, 2005 9:34 am
Location: Calgary, Canada

Post by Darko » Tue Jan 02, 2007 9:18 am

Duh, of course it works :)
(because matrix multiplication is associative)

I still have to deal with that inverse, but I will figure it out eventually.

EonStrife
New poster
Posts: 4
Joined: Thu Jan 04, 2007 3:08 pm

Post by EonStrife » Thu Jan 04, 2007 3:16 pm

The formula I am using is almost similar :
Sum = A * ((A^n)-I) * (A-I)^(-1)

As for the matrix inversion, I have done that (but not sure if It can handle that singular thingy). The problem is, it involves division and floating point, so the final result sometimes off -1 or 1 from desired result (Hello, WA) :(

ibroker
New poster
Posts: 18
Joined: Tue Nov 08, 2005 6:38 pm

11149 Power of Matrix

Post by ibroker » Mon Jan 08, 2007 5:33 am

I use operator overloading calculating matrix.
And I divide equation like this.
A^1+A^2+...A^(2k) = A+A^2+...A^k + A^k(A+A^2..._A^k)
It throws me TLE.
What can I do to cut running time?

Code: Select all


[/code]
Last edited by ibroker on Wed Jan 10, 2007 7:17 am, edited 1 time in total.

Jan
Guru
Posts: 1334
Joined: Wed Jun 22, 2005 10:58 pm
Location: Dhaka, Bangladesh
Contact:

Post by Jan » Tue Jan 09, 2007 12:15 am

Code: Select all

a = calc(k/2);
You are calling it for every recursive call from divide. Try to store the values. Suppose you have stored A^k, then you can find A^(2k) and A^(2k+1) easily and efficiently (Just use the A^k you have stored already). I used this approach and got accepted.
Ami ekhono shopno dekhi...
HomePage

ibroker
New poster
Posts: 18
Joined: Tue Nov 08, 2005 6:38 pm

Re: 11149 TLE..?

Post by ibroker » Wed Jan 10, 2007 7:16 am

Thanks ! I got Accepted

Sanny
Learning poster
Posts: 78
Joined: Sat Feb 14, 2004 3:59 pm
Location: BUET
Contact:

Post by Sanny » Thu Jan 11, 2007 5:27 pm

Are the following i/o correct?

Input:

Code: Select all

1 5
5

5 512
1 2 3 4 5
6 7 8 9 10
11 12 13 14 15
16 17 18 19 20
21 22 23 24 25

4 15
1 2 5 9
8 5 6 7
3 4 8 9
1 3 4 6

3 4095
1 2 3
4 5 6
7 8 9

0 0
Output:

Code: Select all

5

6 2 8 4 0
6 7 8 9 0
6 2 8 4 0
6 7 8 9 0
6 2 8 4 0

8 0 9 4
3 8 9 9
9 5 1 6
0 7 4 5

5 8 1
8 3 8
1 8 5

Thanks in advance.

coolguy
New poster
Posts: 30
Joined: Tue Oct 17, 2006 5:59 pm

Post by coolguy » Thu Jan 11, 2007 7:12 pm

yes they are correct
In good company

Sanny
Learning poster
Posts: 78
Joined: Sat Feb 14, 2004 3:59 pm
Location: BUET
Contact:

Post by Sanny » Thu Jan 11, 2007 9:10 pm

Some more i/o to verify:

Input:

Code: Select all

6 524289
1 2 3 4 5 6
6 5 4 3 2 1
1 3 5 2 4 6
1 3 6 2 5 4
6 2 5 1 3 4
3 4 5 6 2 1 

10 524287
45 46 132 75 45 62 21 19 38 67
45 46 132 75 45 62 21 19 38 67
45 46 132 75 45 62 21 19 38 67
45 46 132 75 45 62 21 19 38 67
45 46 132 75 45 62 21 19 38 67
45 46 132 75 45 62 21 19 38 67
45 46 132 75 45 62 21 19 38 67
45 46 132 75 45 62 21 19 38 67
45 46 132 75 45 62 21 19 38 67
45 46 132 75 45 62 21 19 38 67

2 524287
1 2
3 4
Output:

Code: Select all

6 0 2 1 9 1
1 3 6 4 6 9
6 6 0 7 8 2
7 9 7 0 7 9
2 3 2 7 5 0
2 3 7 5 4 8

5 6 2 5 5 2 1 9 8 7
5 6 2 5 5 2 1 9 8 7
5 6 2 5 5 2 1 9 8 7
5 6 2 5 5 2 1 9 8 7
5 6 2 5 5 2 1 9 8 7
5 6 2 5 5 2 1 9 8 7
5 6 2 5 5 2 1 9 8 7
5 6 2 5 5 2 1 9 8 7
5 6 2 5 5 2 1 9 8 7
5 6 2 5 5 2 1 9 8 7

9 0
5 4


User avatar
sohel
Guru
Posts: 856
Joined: Thu Jan 30, 2003 5:50 am
Location: New York

Post by sohel » Thu Jan 11, 2007 10:17 pm

These are also correct !!

Sanny
Learning poster
Posts: 78
Joined: Sat Feb 14, 2004 3:59 pm
Location: BUET
Contact:

Post by Sanny » Fri Jan 12, 2007 8:04 am

Ok now it is time to post the code i think. I've spent many hours with it. This program gets WA in ~.027 seconds.

Code: Select all

#pragma warning(disable:4786)
#include<vector>
#include<cstdio>
#include<map>
#include<algorithm>
#include<cassert>
using namespace std;

#define rep(i,n) for(i=0;i<(n);i++)

struct matrix
{
	int a[42][42];
}A,I;
int p;
map<int,matrix>mp;//stores different A^n
vector<int>v,vs;
map<int,matrix>S;//S[n] = A + A^2 + ... + A^n

matrix Sum(matrix a,matrix b)
{
	int i,j;
	matrix c;
	rep(i,p) rep(j,p) c.a[i][j] = (a.a[i][j] + b.a[i][j])%10;
	return c;
}

matrix Multiply(matrix a,matrix b)
{
	int i,j,k;
	matrix c;
	rep(i,p) rep(j,p) {
		c.a[i][j] = 0;
		rep(k,p) c.a[i][j] = (c.a[i][j] + a.a[i][k] * b.a[k][j])%10;
	}
	return c;
}

matrix Pow(int n)//finds A^n
{
	if(n==0) return I;
	if(n==1) return A;
	if(mp.find(n) != mp.end()) return mp[n];

	if(n&1) return mp[n] = Multiply(A,Pow(n-1));
	else {
		matrix n2;
		n2 = Pow(n/2);
		return mp[n] = Multiply(n2,n2);
	}
}


//S[n] = S[n/2] + A^(n/2) * S[n/2]; when n is even
//S[n] = A + A * S[n-1]; when n is odd
matrix find_sum(int n)
{
	if(n==0) return I;
	if(n==1) return A;
	if(S.find(n) != S.end()) return S[n];

	if(n&1) return S[n] = Sum(A,Multiply(A,find_sum(n-1)));
	else {
		matrix n2;
		n2 = find_sum(n/2);
		return S[n] = Sum(n2,Multiply(Pow(n/2),n2));
	}
}

//this is used to find which S[n] & A[n] we need in the calculation
void Dummy(int n)
{
	vs.push_back(n);//keeps the S[n]s
	if(n<2) return;
	if(n&1) Dummy(n-1);
	else {
		v.push_back(n/2);//keeps the A[n]s
		Dummy(n/2);
	}
}

int main()
{
	int i,j,n;
	matrix r;
	//p is the size of matrix & n is the power
	while(scanf(" %d %d",&p,&n)==2 && p) {
		assert(p>0 && n>0);
		
		//identity matrix
		rep(i,p) rep(j,p) I.a[i][j] = 0;
		rep(i,p) I.a[i][i] = 1;

		rep(i,p) rep(j,p) {
			scanf("%d",&A.a[i][j]);
			A.a[i][j]%=10;
		}

		mp.clear(); v.clear(); S.clear(); vs.clear();
		mp[0] = I, mp[1] = A;
		S[0] = I, S[1] = A;
		Dummy(n);
		sort(v.begin(),v.end());
		rep(i,v.size()) mp[v[i]] = Pow(v[i]);
		sort(vs.begin(),vs.end());
		rep(i,vs.size()) S[vs[i]] = find_sum(vs[i]);

		r = find_sum(n);
		for(i=0;i<p;i++,printf("\n")) for(printf("%d",r.a[i][0]),j=1;j<p;j++) printf(" %d",r.a[i][j]);
		printf("\n");
	}
	return 0;
}

User avatar
sohel
Guru
Posts: 856
Joined: Thu Jan 30, 2003 5:50 am
Location: New York

Post by sohel » Fri Jan 12, 2007 4:20 pm

I just commented out these five lines

Code: Select all

      Dummy(n); 
      sort(v.begin(),v.end()); 
      rep(i,v.size()) mp[v[i]] = Pow(v[i]); 
      sort(vs.begin(),vs.end()); 
      rep(i,vs.size()) S[vs[i]] = find_sum(vs[i]); 

And it gets AC.

But nothing seems to be wrong with the Dummy()... strange... :-?

Sanny
Learning poster
Posts: 78
Joined: Sat Feb 14, 2004 3:59 pm
Location: BUET
Contact:

Post by Sanny » Fri Jan 12, 2007 4:58 pm

When I remove all the 5 lines, I get AC in 1.473 sec.
If I remove last 2 lines, I get WA in 0.795 sec.
And if I remove 3rd & 4th line only, I get WA in 0.021 sec.

Can anybody explain these?

Post Reply

Return to “Volume 111 (11100-11199)”