Fast power of ACwing scintillation matrix

Posted by vanessa123 on Thu, 13 Jan 2022 20:18:06 +0100

Digression: it was supposed to be posted on Acwing, but I really won't use MARKDOWN on Acwing to post it on CSDN

Matrix fast power

Briefly introduce the following matrix fast exponents (also help yourself to review the following, and those who have already learned to skip by themselves)

Principle: A recursive matrix C is mentioned from the recursive relationship. The initial state is A. each execution is equivalent to (note that the C matrix is in the front)
A = C ∗ A A=C*A A=C∗A
The result matrix B after n times of execution is equal to

B = C n ∗ A B=C^n*A B=Cn∗A

It's natural to think of using fast powers to accelerate this process
The algorithm is divided into two parts: how to get the recursive matrix and the fast power of the matrix, first introduced from the simple fast power

Counting

The principle of fast power will not be repeated. The basic course of y total algorithm is talked about, and you can baidu by yourself. The general fast exponents are as follows:
##Ordinary fast power (c + +)

int ksm(int a,int b)
{
    int res=1;
    while(b)
    {
        if(b&1) 
            res=res*a;
        a=a*a;
        b>>=1;
    }
    return res;
}

What we need to do in the fast power of matrix is to transform the multiplication between integers into the multiplication between matrices

Matrix multiplication is only mentioned here. The number of columns of the first matrix is equal to the number of rows of the second matrix, such as
two × Matrix A of 3 can only be associated with 3 × N is multiplied by matrix B, and the resulting matrix C is 2 × n

There are two ways to change

1. Define the function that implements the matrix function;
2. Names between overloaded matrices;
I use the second one here
##Stick your own board below

struct mat{
	int n, m;
	ll c[20][20];
    /*
        There can be no debug function
        It's just because I don't have enough strength
        Every time I write fast exponents, I have to debug for half a day
    */
	void debug()
	{
		puts("");
		for (int i = 1; i <= n;++i)
		{
			for (int j = 1; j <= m;++j)
				printf("%d ", c[i][j]);
			puts("");
		}
	}
};


mat operator*(mat &x,mat &y)
{
	mat res;
	res.n = x.n;
	res.m = y.m;
	for (int i = 1; i <= res.n;++i)
		for (int j = 1; j <= res.m;++j)
			res.c[i][j] = 0;
	for (int i = 1; i <= res.n;++i)
	{
		for (int j = 1; j <= res.m;++j)
		{
			for (int k = 1; k <= x.m;++k)
			{
				res.c[i][j] = (res.c[i][j] + x.c[i][k] * y.c[k][j] % mod) % mod;
			}
		}
	}
	return res;
}

In this way, you can directly use the board of fast power above to realize the fast power between matrices. The difference is that the data type has changed

The recursive matrix is obtained

This is the difficulty in the algorithm. Through the matrix fast power template above, as long as you get the correct recursive matrix, just A.
Let's start with a simple one: fiboracci sequence, recursive formula:
{ f ( n ) = f ( n − 1 ) + f ( n − 2 ) , n > = 2 f ( 1 ) = 1 , f ( 2 ) = 1 \begin{cases} f(n)=f(n-1)+f(n-2),n>=2\\f(1)=1,f(2)=1 \end{cases} {f(n)=f(n−1)+f(n−2),n>=2f(1)=1,f(2)=1​
The first two terms are both 1, so the initial matrix A:
( 1 1 ) \begin{pmatrix} 1\\1 \end{pmatrix} (11​)
According to the recurrence relationship, the current term n is the sum of the first two terms n-1 and n-2, and the N-1 term remains unchanged. Therefore, the recurrence matrix C is:
( 1 1 0 1 ) \begin{pmatrix} 1 &1\\ 0&1 \end{pmatrix} (10​11​)
It may be difficult for beginners to understand here. You can try to execute C*A once, that is
( 1 ∗ 1 + 1 ∗ 1 = 2 0 ∗ 1 + 1 ∗ 1 = 1 ) \begin{pmatrix} 1*1+1*1=2\\0*1+1*1=1 \end{pmatrix} (1∗1+1∗1=20∗1+1∗1=1​)
It conforms to the recursive formula.
The above is a brief introduction to the fast power of matrix.
Back to the topic

Recursive matrix

The following conclusions can be drawn from the title:

In the m second, the value of the nth bulb is equal to the m-1 second, and the value of the n-1 bulb is different from the value of the m-1 second and the nth bulb, or
Equivalence and
In the m second, the value of the nth bulb is equal to the m-1 second, and the addition of the value of the n-1 bulb and the value of the m-1 second and the nth bulb in the sense of touch 2
Write recursively:
f ( m , n ) = f ( m − 1 , n ) ⊕ f ( m − 1 , n − 1 ) ⇒ f ( m , n ) = f ( m − 1 , n ) + f ( m − 1 , n − 1 ) ( m o d 2 ) f(m,n)=f(m-1,n) \oplus f(m-1,n-1)\Rightarrow \\ f(m,n)=f(m-1,n)+f(m-1,n-1)\quad(mod\quad2) f(m,n)=f(m−1,n)⊕f(m−1,n−1)⇒f(m,n)=f(m−1,n)+f(m−1,n−1)(mod2)

The recursive matrix is obtained
( 1 0 0 . . . 1 1 1 0 . . . 0 0 1 1 . . . 0 ⋮ ⋮ ⋮ ⋮ ⋮ 0 0 . . . 1 1 ) \begin{pmatrix} 1&0&0&...&1\\ 1&1&0&...&0\\ 0&1&1&...&0\\ \vdots&\vdots&\vdots&\vdots&\vdots\\ 0&0&...&1&1 \end{pmatrix} ⎝⎜⎜⎜⎜⎜⎛​110⋮0​011⋮0​001⋮...​.........⋮1​100⋮1​⎠⎟⎟⎟⎟⎟⎞​
The initial matrix is given by the title

Code entry

#include<iostream>
#include<cstring>
using namespace std;
typedef long long ll;
const int N=20;

int n;
ll k;

struct mat{
    int n,m;
    int c[N][N];

    void debug()
    {
        for(int i=1;i<=n;++i)
        {
            for(int j=1;j<=m;++j)
                printf("%d ",c[i][j]);
            puts("");
        }
    }

};

mat operator*(mat &x,mat &y)
{
    mat res;
    res.n=x.n;
    res.m=y.m;
    for(int i=1;i<=res.n;++i)
        for(int j=1;j<=res.m;++j)
            res.c[i][j]=0;
    for(int i=1;i<=res.n;++i)
        for(int j=1;j<=res.m;++j)
            for(int k=1;k<=x.m;++k)
                res.c[i][j]=(res.c[i][j]^(x.c[i][k]*y.c[k][j]));
    return res;
}

void ksm()
{
    cin>>n>>k;
    mat res,a;
    res.n=n;
    res.m=1;
    a.n=a.m=n;
    for(int i=1;i<=n;++i)
        for(int j=1;j<=n;++j)
            a.c[i][j]=0;
    a.c[1][1]=1;
    a.c[1][n]=1;
    for(int i=2;i<=n;++i)
    {
        a.c[i][i-1]=1;
        a.c[i][i]=1;
    }
    for(int i=1;i<=n;++i)
        cin>>res.c[i][1];
    while(k)
    {
        if(k&1)
            res=a*res;
        a=a*a;
        k>>=1;
    }
    res.debug();
}

int main()
{
    ksm();
    return 0;
}

Later: you can directly use XOR operation to overload the multiplication between matrices (bit operation should be faster, but the amount of data is too small to reflect)
It's 3:06 on January 14, 2022. Come on, everyone.

Topics: C++ linear algebra