從0開始學算法--數學(4.1矩陣)


1,概念

以下詞條解釋來自百度百科:代數,代數系統,線性代數,矩陣

代數

  代數是研究數、數量、關系、結構與代數方程(組)的通用解法及其性質的數學分支。初等代數一般在中學時講授,介紹代數的基本思想:研究當我們對數字作加法或乘法時會發生什么,以及了解變量的概念和如何建立多項式並找出它們的根。代數的研究對象不僅是數字,而是各種抽象化的結構。在其中我們只關心各種關系及其性質,而對於“數本身是什么”這樣的問題並不關心。常見的代數結構類型有群、、域、模、線性空間等。

代數系統

  非空集合A和A上k個一元或二元運算f1,f2,…,fk組成的系統稱為一個代數系統,簡稱代數,記作(A,f1,f2,…,fk)。由定義可知,一個代數系統需要滿足下面3個條件:(1)有一個非空集合A;(2)有一些建立在集合A上的運算;(3)這些運算在集合A上是封閉的。有的書上對代數系統定義時不要求運算的封閉性,而是把具有封閉性的代數系統定義為一個新的概念-廣群

線性代數

  線性代數是數學的一個分支,它的研究對象是向量向量空間(或稱線性空間),線性變換和有限維的線性方程組。向量空間是現代數學的一個重要課題;因而,線性代數被廣泛地應用於抽象代數泛函分析中;通過解析幾何,線性代數得以被具體表示。線性代數的理論已被泛化為算子理論。由於科學研究中的線性模型通常可以被近似為線性模型,使得線性代數被廣泛地應用於自然科學和社會科學中。

矩陣

  在數學中,矩陣(Matrix)是一個按照長方陣列排列的復數實數集合 [1]  ,最早來自於方程組系數常數所構成的方陣

2,矩陣的基本運算

矩陣類

矩陣的信息有行數列數和矩陣的內容

struct matrix{
    int n,m;
    long long a[N][M];
    matrix(){//    初始化2*2的單位矩陣
        n=m=2;
        for(int i=0;i<n;i++){
            for(int j=0;j<m;j++){
                a[i][j]=0;
            }
        }
        for(int i=0;i<n;i++){
            a[i][i]=1;
        }
    }
    void clear(){
        //memset(a,0,sizeof(a));
        for(int i=0;i<n;i++){
            for(int j=0;j<m;j++){
                a[i][j]=0;
            }
        }
    }
}

加和減

矩陣的加減運算要抱枕兩個矩陣的行數和列數相同

matrix operator+(const matrix &b)const{
    matrix tmp;
    tmp.clear();
    for(int i=0;i<n;i++){
        for(int j=0;j<m;j++){
            tmp.a[i][j]=a[i][j]+b.a[i][j];
        }
    }
    return tmp;
}

乘法

設矩陣乘法的到的結果是矩陣S,則矩陣S的第i行第j個元素被乘數的第i行元素乘以乘數的第j列元素的求和得到的,也就是說被乘數的列數乘數的行數相等

matrix operator*(const matrix &b)const{
    matrix tmp;
    tmp.clear();
    for(int i=0;i<n;i++){
        for(int j=0;j<b.m;j++){
            for(int k=0;k<m;k++){
                tmp.a[i][j]=((a[i][k]*b.a[k][j])%mod+tmp.a[i][j])%mod;
            }
        }
    }
    return tmp;
}

因為乘法對矩陣有要求,所以矩陣乘法並不滿足交換律,只滿足結合律

除法

矩陣並沒有除法運算,但是可以求逆,類似在初等數學中數的乘法逆元是他的倒數。矩陣的逆也滿足這樣的條件,1.矩陣S的逆的逆是矩陣S本身,2.除以矩陣S等於乘以矩陣S的逆

矩陣求逆代碼:

vector<double> operator * (vector<double> a,double b){
    int n=a.size();
    vector<double> res(n,0);
    for(int i=0;i<n;i++){
        res[i]=a[i]*b;
    }
    return res;
}

vector<double> operator - (vector<double> a,vector<double> b){
    int n=a.size();
    vector<double> res(n,0);
    for(int i=0;i<n;i++){
        res[i]=a[i]-b[i];
    }
    return res;
}

void get_ni(vector<double>a[],vector<double>b[],int n){//b為初等變換所需要的單位矩陣
    for(int i=0;i<n;i++){
        b[i]=vector<double>(n,0);
    }
    for(int i=0;i<n;i++){//初始化單位矩陣
        b[i][i]=1;
    }
    for(int i=0;i<n;i++){//初等變換
        for(int j=i;j<n;j++){
            if(fabs(a[j][i])>0){
                swap(a[i],a[j]);
                swap(b[i],b[j]);
                break;
            }
        }
        b[i]=b[i]*(1.0/a[i][i]);
        a[i]=a[i]*(1.0/a[i][i]);
        for(int j=0;j<n;j++){
            if(j!=i&&fabs(a[j][i]>0)){
                b[j]=b[j]-b[i]*a[j][i];
                a[j]=a[j]-a[i]*a[j][i];
            }
        }i
    }
}

3.高斯消元

略。。。

4.常系數線性齊次遞推

在求解斐波那契數列問題中,如果使用記憶化的技巧,在線性時間內就可以求解第n項斐波那契數是多少。

在矩陣中這個問題如何解決??

如果可以構造如圖所示的矩陣,那么連續給矩陣乘以n個這樣的矩陣就可以得到fn

 

 又因為矩陣滿足結合律,所以可以用快速冪的方式,除去矩陣乘法的時間復雜度,O(logn)就可以得到斐波那契數列的第n項。

對於常系數線性齊次遞推得到如下結論:

 

 代碼:

const int N=10,M=10;
long long a[10];
long long f[10];

struct matrix{
    int n,m;
    long long a[N][M];
    matrix(){//    初始化2*2的單位矩陣
        n=m=2;
        for(int i=0;i<n;i++){
            for(int j=0;j<m;j++){
                a[i][j]=0;
            }
        }
        for(int i=0;i<n;i++){
            a[i][i]=1;
        }
    }
    void clear(){
        //memset(a,0,sizeof(a));
        for(int i=0;i<n;i++){
            for(int j=0;j<m;j++){
                a[i][j]=0;
            }
        }
    }
    matrix operator+(const matrix &b)const{
        matrix tmp;
        tmp.clear();
        for(int i=0;i<n;i++){
            for(int j=0;j<m;j++){
                tmp.a[i][j]=a[i][j]+b.a[i][j];
            }
        }
        return tmp;
    }
    matrix operator*(const matrix &b)const{
        matrix tmp;
        tmp.clear();
        for(int i=0;i<n;i++){
            for(int j=0;j<b.m;j++){
                for(int k=0;k<m;k++){
                    tmp.a[i][j]=((a[i][k]*b.a[k][j])+tmp.a[i][j]);
                }
            }
        }
        return tmp;
    }
};

matrix pow2(matrix A,int n){
    matrix B;
    while(n>0){
        if(n&1)B=B*A;
        A=A*A;
        n/=2;
    }
    return B;
}

long long solve(long long a[],long long f[],int n,int m){
    matrix A,B;
    A.clear();
    B.clear();
    A.n=B.n=n;
    A.m=n;
    B.m=1;
    for(int i=0;i<n-1;i++){
        A.a[i+1][i]=1;
    }
    for(int i=0;i<n;i++){
        B.a[i][0]=f[n-i-1];
        A.a[0][i]=a[n-i-1];
    }
    A=pow2(A,m);
    B=A*B;
    return B.a[n-1][0];
}

int main(){//求斐波那契數列第m項
    a[0]=1,a[1]=1;      //系數
    f[0]=0,f[1]=1;      //前n項
    int m;
    scanf("%d",&m);
    printf("%lld\n",solve(a,f,2,m));
    return 0;
}

一點心得:(走在路上啃饅頭想到的(滑稽)

如果我們的遞推式不僅僅是fn=an-1*fn-1+an-2*fn-2+………+a0*f0;這樣簡單的格式.而是fn=an-1*fn-1+an-2*fn-2+………+a0*f0+c;c是一個常數。那么可以構造如下的矩陣解決一個問題。

 

構造矩陣的方法或許沒有多么巧妙,我卻收獲頗多。

1.矩陣乘法本身其實就在闡述一定的遞推關系。

2.如果遇到下圖所示矩陣,可分解為兩個矩陣。

 

3.如果可以構造出f(n)=n的矩陣,則可以輕易得出fn=an-1*fn-1+an-2*fn-2+………+a0*f0+n的矩陣遞推式

同理如果可以構造f(n)=an的矩陣,就可以輕易得出fn=an-1*fn-1+an-2*fn-2+………+a0*f0+an的矩陣遞推式

結論:構造出的任意矩陣可以嵌套出一個新的遞推式。

 例:

5.遞推(第二彈)

例題:POJ3734 Blocks(自我感覺本題闡述了矩陣本身就代表着一定的遞推關系

題意:給出n個排成一列的方塊,用紅、藍、綠、黃四種顏色給它們染色,要求出染成紅色的方塊個數和染成綠色的方塊個數同時為偶數的染色方案的個數模10007的值為多少。

從左邊開始考慮,染第i個格子的時候有三種情況:

1. 之前i-1個方塊中,紅色數和綠色數都是偶數,想要滿足條件,就要染藍色或黃色。

2. 之前i-1個方塊中,紅色數和綠數中有一個為奇數,想要滿足條件就要染奇數個對應的顏色。

3. 之前i-1個方塊中,紅色的方塊個數和綠色的方塊個數均不是偶數,這時無法滿足條件。

如果用a[],b[],c[]三個數組表示三種狀態。

那么可以得到三個遞推式

ai=2*ai-1+bi-1

bi=2*ai-1+2*bi-1+2*ci-1

ci=bi-1+2*ci-1

那么可以得到如下矩陣:

 

 矩陣快速冪求解即可得到答案.

代碼:

#include <iostream>
#include <algorithm>
#include <cstdio>
#include <cstring>

using namespace std;

const int N=10,M=10;
long long a[10];
long long f[10];
const int mod=10007;


struct matrix{
    int n,m;
    long long a[N][M];
    matrix(){//    初始化2*2的單位矩陣
        n=m=3;
        for(int i=0;i<n;i++){
            for(int j=0;j<m;j++){
                a[i][j]=0;
            }
        }
        for(int i=0;i<n;i++){
            a[i][i]=1;
        }
    }
    void clear(){
        //memset(a,0,sizeof(a));
        for(int i=0;i<n;i++){
            for(int j=0;j<m;j++){
                a[i][j]=0;
            }
        }
    }
    matrix operator+(const matrix &b)const{
        matrix tmp;
        tmp.clear();
        for(int i=0;i<n;i++){
            for(int j=0;j<m;j++){
                tmp.a[i][j]=a[i][j]+b.a[i][j];
            }
        }
        return tmp;
    }
    matrix operator*(const matrix &b)const{
        matrix tmp;
        tmp.clear();
        for(int i=0;i<n;i++){
            for(int j=0;j<b.m;j++){
                for(int k=0;k<m;k++){
                    tmp.a[i][j]=((a[i][k]*b.a[k][j])%mod+tmp.a[i][j])%mod;
                }
            }
        }
        return tmp;
    }
};

matrix pow2(matrix A,int n){
    matrix B;
    while(n>0){
        if(n&1)B=B*A;
        A=A*A;
        n/=2;
    }
    return B;
}


long long solve(int n,int m){
    matrix A,B;
    A.clear();
    B.clear();
    A.n=B.n=n;
    A.m=n;
    B.m=1;
    A.a[0][0]=2,A.a[0][1]=1,A.a[0][2]=0;
    A.a[1][0]=A.a[1][1]=A.a[1][2]=2;
    A.a[2][0]=0,A.a[2][1]=1,A.a[2][2]=2;
    B.a[0][0]=1,B.a[0][1]=0,B.a[0][2]=0;
    A=pow2(A,m);
    B=A*B;
    return B.a[0][0];
}

int main(){

    int n;
    scanf("%d",&n);
    for(int i=0;i<n;i++){
        int m;
        scanf("%d",&m);
        printf("%d\n",solve(3,m));
    }
    return 0;
}

 

6.十進制快速冪

例:2019牛客暑期多校訓練營(第五場)B題

 

 

 

 

 

 

 

 

 代碼:

#include<algorithm>
#include<iostream>
#include<cstdlib>
#include<fstream>
#include<cstring>
#include<bitset>
#include<cstdio>
#include<time.h>
#include<deque>
#include<queue>
#include<stack>
#include<cmath>
#include<map>
#include<set>
 
using namespace std;
 
long long mod;
const int N=2,M=2;
const int maxnn=1e6+10;
long long a[10];
long long b[10];
char t[maxnn];
 
struct matrix{
    long long a[N][M];
    matrix(){
        for(int i=0;i<2;i++){
            for(int j=0;j<2;j++){
                a[i][j]=0;
            }
        }
        for(int i=0;i<2;i++){
            a[i][i]=1;
        }
    }
    void clear(){
        for(int i=0;i<2;i++){
            for(int j=0;j<2;j++){
                a[i][j]=0;
            }
        }
    }
    matrix operator+(const matrix &b)const{
        matrix tmp;
        tmp.clear();
        for(int i=0;i<2;i++){
            for(int j=0;j<2;j++){
                tmp.a[i][j]=a[i][j]+b.a[i][j];
            }
        }
        return tmp;
    }
    matrix operator*(const matrix &b)const{
        matrix tmp;
        tmp.clear();
        for(int i=0;i<2;i++){
            for(int j=0;j<2;j++){
                for(int k=0;k<2;k++){
                    tmp.a[i][j]=((a[i][k]*b.a[k][j])%mod+tmp.a[i][j])%mod;
                }
            }
        }
        return tmp;
    }
};
 
matrix pow2(matrix A,int n){
    matrix B;
    while(n>0){
        if(n&1)B=B*A;
        A=A*A;
        n/=2;
    }
    return B;
}
 
int main(){
    for(int i=0;i<2;i++){
        scanf("%lld",&a[i]);
    }
    for(int i=1;i>=0;i--){
        scanf("%lld",&b[i]);
    }
    getchar();
    scanf("%s",t);
    getchar();
    scanf("%lld",&mod);
    matrix A;
    A.clear();
    A.a[1][0]=1;
    for(int i=0;i<2;i++){
        A.a[0][i]=b[2-i-1];
    }
    int len=strlen(t);
    matrix D;
    for(int i=len-1;i>=0;i--){
        int x=(t[i]-'0');
        D=D*pow2(A,x);
        A=pow2(A,10);
    }
    long long num=((D.a[1][0]*a[1])%mod+(D.a[1][1]*a[0])%mod)%mod;
    printf("%lld\n",num);
    return 0;
}

 7.補:

一:poj-3233

給定n*n的矩陣A和正整數k和m,求矩陣A的冪的和。S=A+A2+A3+…………+Ak

輸出S的各元素對mod取余后的答案。

n<30,k<109,m<104.

思路一:按線性求出所有A的x次方,加起來,時間復雜度為O(n3k)

思路二:

矩陣的也可以看成一個元素,可以把矩陣元素放在矩陣中加速這個過程。

令Sk=I+A+……+Ak-1(I是單位矩陣)

 

 

這樣就可以在O(n3logk)的時間復雜度內解決問題。

 


免責聲明!

本站轉載的文章為個人學習借鑒使用,本站對版權不負任何法律責任。如果侵犯了您的隱私權益,請聯系本站郵箱yoyou2525@163.com刪除。



 
粵ICP備18138465號   © 2018-2026 CODEPRJ.COM