stick_矩阵乘法

好吧， 重新写一写矩阵乘法。

考虑斐波拉契数列 Fi = Fi-1+ Fi-2 

那么考虑一个1*2的矩阵{ Fi,   Fi-1 } 乘以一个 2*2的矩阵， 成为{ Fi+1 ,  Fi}

考虑矩阵乘法的定义， 这个一维数组的第一个数即为 原一位数组的每一个数 乘以 2*2矩阵的第一列数；

第二个数即 原一位数组的每一个数 乘以 2*2矩阵的第二列数；

以此类推。

所以因为Fi+1 = Fi + Fi-1 所以 第一列数为 {1, 1}

Fi= Fi  所以第二列数为{ 0, 1}

我们可以通过递推式来确定矩阵。

再考虑另一个操作。

A=  A+B-C

B= 2A+C

C= B-A

把这三个数进行n次操作。

考虑一维矩阵(A, B, C)- 乘以一个矩阵， 使之成为新的A,B,C

这个矩阵为

{ 1， 2， -1，

  1，  0， 1

 -1， 1 ， 0}

这是根据系数得到的。

再看看一道题。Fibonacci的变种。

Fi= Fi-1 + Fi-2 + i

首先考虑状态有几个数~

两个数是不行的~~~~那我们设状态为 {Fi , Fi-1, i} 能不能推出 {Fi+1, Fi, i+1}

这里包括了两个操作： Fi= Fi-1 + Fi-2 + i;    i+1= i + 1;那我们发现后一个是求不出来的~~因为不能实现自增操作~~那么只能在加多一个数“1”。

所以状态变成了{Fi, Fi-1, i, 1} 

于是， 矩阵就是：

{

1， 1， 0， 0

1， 0， 0， 0

1， 0， 1， 0

0， 0， 1， 1

}

言归正传， 看看题目。 时限5s。<最好是1s过>

【问题描述】

CD拿出来了好多1×1×2的积木，准备按像这样的俯视图的方式堆积木：

XXX

XOX

XXX {O是空的}

正如你看到的，是一个3×3中间是空的正方形！

例如，对于第一层，其中一种方式可以是这样的：

AAB

COB

CDD

现在，CD计划堆一个高度为N层的积木（除了中心外，其他地方不能留空），那么，他又多少种堆积木的方法呢？

【输入】

输入一行，包含一个整数N(1<=1000000000)。

【输出】

输出所有情况的数量对1000000007取模后的结果。

 

这道题递推是比较明显的递推， 设F[i, j]为到了第i层, i层以下以被填满，第i层 还没有被填满， 但已经被i-1层竖起来的积木占用了状态为j的位置。

由于只有8个位置，总状态数为2^8=256个。

递推式可以枚举每一种状态, 再用深搜摆满, 记录每一种合法的状态的竖着摆放的状态。即 x 可以 推向 y。

则 F[i+1, k] += F[ i, j ]  (状态j可以推向状态k)

由于i维巨大， 考虑矩阵乘法。

将F[i] 【一个一维矩阵】  乘以一个矩阵【2维】 得到 F[i+1]

矩阵可以由递推式得到。 自行YY。

由于矩阵乘法有结合律。所以可以先计算这个矩阵 G 【即上面求出的2维矩阵】 的 n次方 算出来， 再与 一维矩阵相乘。

于是可以使用快速幂。

x^y = (x^ y/2)^2 (y mod 2=0)

x^y = (x^ y/2)^2 * x (y mod 2=1)【此处为整除】

矩阵乘法的时间复杂度为N^3

于是时间为 O(256^3 * log N )

The End.

超级优化：

把无用的状态去掉。【把从0开始搜索到的状态才记录下来】

于是时间可以优化到(70^3 * log N)

再贴一下恶心的代码~~~

#include<cstdio>
#include<cstring>
using namespace std;

#define N 256
#define Mod 1000000007

const int pair[8][2]={1,3, 0, 2, 1, 4, 0, 5, 2, 7, 3, 6, 5, 7, 4, 6};

typedef long long int64;
typedef int64 mat[N][N];

mat yuan, now, tp;
int n, O;
int lab[N], que[N];

bool f[8], g[8], w[N];

void cha(int x){
     memset(f, 0, sizeof(f));
     for(int i=0; i<=7; i++, x>>=1) f[i]= (x & 1);
     }
     
int coding(bool f[]){
    int ret=0;
    for(int i=7; i>=0; i--) {ret=ret*2; if(f[i]) ret++;}
    return ret;
}

void dfs(int x, int y){
     if(y==8){
             int code= coding(g);
             yuan[lab[code]][x]++;
             //yuan[x][code]++;
             return;
             }
     if(f[y]){ dfs(x, y+1); return; }
     int p= pair[y][0], q=pair[y][1], ty=y+1;
     if(!f[p]){   f[y]=f[p]=1;   dfs(x, y+1);   f[y]=f[p]=0;}
     if(!f[q]){   f[y]=f[q]=1;   dfs(x, y+1);   f[y]=f[q]=0;}
     f[y]=g[y]=1;   dfs(x, y+1);   f[y]=g[y]=0;
     return;
     }

void expand(int y){ // 我竟然写了和一个dfs一样的过程来求有用的状态！！！
     if(y==8){
             int code= coding(g);
             if(lab[code]<0) {que[++O]=code; lab[code]=O;}
             return;
             }
     if(f[y]){ expand(y+1); return; }
     int p= pair[y][0], q=pair[y][1], ty=y+1;
     if(!f[p]){   f[y]=f[p]=1;   expand(y+1);   f[y]=f[p]=0;}
     if(!f[q]){   f[y]=f[q]=1;   expand(y+1);   f[y]=f[q]=0;}
     f[y]=g[y]=1;   expand(y+1);   f[y]=g[y]=0;
     return;
     }

void graph(){
     memset(lab, 255, sizeof(lab));
     O=0;
     que[0]=0;
     lab[0]=0;
     int t=0;
     while(t<=O){
                      cha(que[t]);
                      memset(g, 0 ,sizeof(g));
                      expand(0);
                      ++t;
                      }
     
     for(int i=0; i<N; i++)if(lab[i]>=0){
            cha(i);
            memset(g, 0 ,sizeof(g));
            dfs(lab[i], 0);
            }
     }
     
void mult(mat a, mat b){
    memset(tp, 0, sizeof(tp));
    for(int i=0;i<=O;i++)
            for(int j=0;j<=O;j++)
                    for(int k=0;k<=O;k++) tp[i][j]= (tp[i][j] + b[i][k] * a[k][j] )% Mod;
    memcpy(now, tp, sizeof(tp));
}
     
void solve(){
     int tn=n, tp=0;
     while(tn>0){ w[++tp]= (tn & 1); tn>>=1; }
     memcpy(now, yuan, sizeof(now));
     for(int i=tp-1; i>0; i--){
             mult(now, now);
             if(w[i]) mult(yuan ,now);
             }
     printf("%d\n", now[0][0]);
     }
     
int main(){
    freopen("stick.in","r",stdin);
    freopen("stick.out","w",stdout);
    scanf("%d",&n);
    graph();
    solve();
    fclose(stdin);fclose(stdout);
    return 0;
}