垒骰子与矩阵快速幂

转自 http://blog.csdn.net/lonverce/article/details/45169285

垒骰子

 

赌圣atm晚年迷恋上了垒骰子，就是把骰子一个垒在另一个上边，不能歪歪扭扭，要垒成方柱体。

经过长期观察，atm 发现了稳定骰子的奥秘：有些数字的面贴着会互相排斥！

我们先来规范一下骰子：1 的对面是 4，2 的对面是 5，3 的对面是 6。

假设有 m 组互斥现象，每组中的那两个数字的面紧贴在一起，骰子就不能稳定的垒起来。 

atm想计算一下有多少种不同的可能的垒骰子方式。

两种垒骰子方式相同，当且仅当这两种方式中对应高度的骰子的对应数字的朝向都相同。

由于方案数可能过多，请输出模 10^9 + 7 的结果。

 

不要小看了 atm 的骰子数量哦～

 

解法一：利用动态规划，速度还是不够快

#include<iostream>
using namespace std;

//compact[i][j]=false代表点数i,j的面存在冲突
bool Compact[7][7];

//Parner[i]=j代表点数为i的面的对立面点数为j
const int Parner[7]={0,4,5,6,1,2,3};
const long long MOD=1000000007;

int main()
{
    long long N;
    int M;
    int s1,s2;
    cin>>N>>M;
    for(int i=0;i<7;i++)
        for(int j=0;j<7;j++)
        Compact[i][j]=true;
    for(int i=0;i<M;i++)
    {
        cin>>s1>>s2;
        Compact[s1][s2]=Compact[s2][s1]=false;
    }
    long long dp[2][7]; //滚动数组
    long long C=4;
    int e=0;
    for(int i=1;i<7;i++)
        dp[e][i]=1;
    //dp[i][j]代表高度为i的,顶面点数为j的叠骰子方案数 
    //在这里忽略每个骰子可以四面转向的情况, 把该情况留到最后乘上去就可以了
    int j,k;
    for(long long i=2;i<=N;i++)
    {
        e=1-e;  //滚动处理
        C=(C*4)%MOD;
        for(j=1;i<7;j++)
        {
            dp[e][j]=0;
            for(k=1;k<7;k++)
                if(Compact[Parner[j]][k])
                    dp[e][j]+=dp[1-e][k];
            dp[e][j]%=MOD;
        }
    }
    int sum=0;
    for(int i=1;i<7;i++)
        sum=(sum+dp[e][i])%MOD;
    sum=(sum*C)%MOD;
    cout<<sum;
    return 0;
}

解法二：矩阵快速幂

传说中用矩阵快速幂做的题。
首先看动态规划的思路。
Dp[i][j]表示高度为i，顶面点数为j的方案数，那么Dp[i][j]就等于i-1高度是所有与j的反面无冲突的方案累加，最后的总方案数还要乘以4^i
在本题的动态规划解法中，以dp[i][j]来表示在高度为i的情况下，骰子柱顶面点数为j的总方案数，然后。dp[i][j]的值将等于前一高度所有方案的选择性累加。选择性累加就类似包含0,1，的矩阵。


接下来，我们需要沿用在DP解法中提到的冲突数组，因为冲突数组实际就是我们的选择性累加，在这里我们把它变成冲突矩阵，并且更正对他的含义

#include<iostream>
#include<vector>
using namespace std;

typedef unsigned int uint;
typedef uint SizeType;  //表示大小
typedef uint IndexType;  //表示下标
typedef long long ValueType;  //表示值
typedef vector<ValueType>RowType;  //表示矩阵的一行
typedef vector<RowType>ContainerType;  //矩阵容器
#define MOD 1000000007

//矩阵类的实现
class Matrix{
    public:
        Matrix( SizeType Row, SizeType Col, ValueType Init_Val=0 ):m_row(Row),m_col(Col){
            m_mtx = new ContainerType( m_row, RowType(m_col,Init_Val) );
        }
        Matrix( const Matrix& Copyright ):m_row(Copyright.m_row),m_col(Copyright.m_col){
            m_mtx = new ContainerType( *Copyright.m_mtx );
        }
        inline SizeType GetRow()const{ return m_row; }
        inline SizeType GetCol()const{ return m_col; }

        // ...矩阵赋值
        const Matrix& operator=( const Matrix& right ){
            *m_mtx = *right.m_mtx;
            m_row = right.m_row;
            m_col = right.m_col;
        }

        // ...获取矩阵中第rI行cI列的元素
        inline const ValueType& operator()( IndexType rI,IndexType cI )const{
            return (*m_mtx)[rI][cI];
        }
        inline ValueType& operator()( IndexType rI, IndexType cI){
            return const_cast<ValueType&>( static_cast<const Matrix&>(*this)(rI,cI) );
        }

        // ...获得一个n阶的单位矩阵
        static Matrix GetIdentityMatrix( SizeType nDimension )
        {
            Matrix I( nDimension,nDimension,0);
            for( IndexType i = 0; i < nDimension; ++i)
                I(i,i) = 1;
            return I;
        }

        // ...矩阵乘法
        const Matrix operator*( const Matrix& right )const{
            Matrix Res( m_row, right.m_col,0);
            if( m_col == right.m_row ){
                const Matrix& left = *this;
                for( IndexType i = 0; i < m_row; ++i)
                    for( IndexType j = 0; j < right.m_col; ++j)
                        for( IndexType t = 0; t < m_col; ++t )
                            Res(i,j) = (Res(i,j)+left(i,t)*right(t,j)%MOD)%MOD;
            }
            return Res;
        }
        const Matrix& operator*=( const Matrix& right ){return *this = (*this)*right;}

        // ...矩阵快速幂
        const Matrix operator^( uint N ) const{
            Matrix Res = GetIdentityMatrix( m_row ),t = *this;
            while( N ){
                if( N&1 ) Res*=t;
                t*=t;
                N>>=1;
            }
            return Res;
        }
        const Matrix& operator^=(uint N){ return *this = (*this)^N; }
    protected:
        SizeType m_row, m_col;  // ...行数, 列数
        ContainerType *m_mtx;   // ...容器
};
//常数快速幂
ValueType Pow(ValueType a,uint N)
{
    ValueType Res=1;
    while(N)
    {
        if(N&1)
            Res=Res*a%MOD;
        a=a*a%MOD;
        N>>=1;
    }
    return Res;
}
int main(int argc, char** argv) {
    Matrix Complict(6,6,1); // ...冲突组合记录
    Matrix Count( 1,6,1 );// ...总数记录
    const IndexType Parner[6] = {3,4,5,0,1,2}; // ...反面记录
    uint N, M, s1, s2;
    cin >> N >> M;
    for( uint i = 0; i < M; ++i){
        cin >> s1 >> s2;
        Complict( Parner[s1],s2 ) = Complict( Parner[s2],s1 ) = 0;
    }
    Complict^=(N-1);
    Count*=Complict;
    ValueType sum = 0;
    for( IndexType i = 0; i < Count.GetCol(); ++i)
        sum = (sum+Count(0,i))%MOD;
    cout << sum*Pow(4,N)%MOD;
    return 0;
}