In mathematics you don't understand things, you just get used to them.

标签 matrices 下的文章

Description

Link.

$$ \sum\prod_{i=1}^{m}F_{a_{i}},(m>0,a_{1},\cdots a_{m}>0,\sum a_{i}=n) $$

Solution

这是一篇不用 $\mathbf{OGF}$ 的题解。

设 $f_{i}$ 为 $i$ 的 $\operatorname{lqp}$ 拆分值。

然后有显然的过不了递推式:

$$ f_{n}=\begin{cases} 1,n=0 \\ \displaystyle \sum_{i=0}^{n-1}F_{n-i}\times f_{i},n\neq0 \end{cases} $$

然后传统艺能错位相减操作一下:

$$ \begin{aligned} f_{n}&=\sum_{i=0}^{n-1}F_{n-i}\times f_{i} \\ f_{n-1}&=\sum_{i=0}^{n-2}F_{n-i-1}\times f_{i} \\ f_{n-2}&=\sum_{i=0}^{n-3}F_{n-i-2}\times f_{i} \end{aligned} \Longrightarrow \begin{aligned} f_{n}-f_{n-1}-f_{n-2}&=\sum_{i=0}^{n-1}F_{n-i}\times f_{i}-\sum_{i=0}^{n-2}F_{n-i-1}\times f_{i}-\sum_{i=0}^{n-3}F_{n-i-2}\times f_{i} \\ f_{n}-f_{n-1}-f_{n-2}&=(F_{2}-F_{1})\times f_{n-2}+F_{1}\times f_{n-1} \end{aligned} \\ \downarrow \\ f_{n}=2f_{n-1}+f_{n-2} $$

递推公式有了,然后矩阵快速幂:

$$ \begin{bmatrix} f_{n} \\ f_{n-1} \end{bmatrix} =\begin{bmatrix} 2f_{n-1}+f_{n-2} \\ f_{n-1} \end{bmatrix} =\begin{bmatrix} f_{n-1} \\ f_{n-2} \end{bmatrix} \times\begin{bmatrix} 2 & 1 \\ 1 & 0 \end{bmatrix} $$

这样就可以做了(吗?):

(code?)

#include <cstdio>
#include <iostream>
#include <cstring>
#include <queue>
#define mod ( 1000000007 )

using namespace std;
typedef long long LL;

template<typename _T, typename _P>
_T qkpow( _T bas, _T one, _P fur ){
    _T res = one;
    while( fur != 0 ){
        if( fur % 2 == ( _P )1 )    res = bas * res;
        bas = bas * bas;
        fur /= 2;
    }
    return res;
}

template<typename _T>
_T add( _T x, _T y ){ if( y >= mod )    y %= mod; x += y; if( x >= mod )    x -= mod; return x; }

struct bigInt : vector<int>{
    bigInt &check( ){
        while( ! empty( ) && ! back( ) ) pop_back( );
        if( empty( ) )    return *this;
        for( unsigned i = 1; i < size( ); ++ i ){ ( *this )[i] += ( *this )[i - 1] / 10; ( *this )[i - 1] %= 10; }
        while( back( ) >= 10 ){ push_back( back( ) / 10 ); ( *this )[size( ) - 2] %= 10; }
        return *this;
    }
    bigInt( int tpN = 0 ){ push_back( tpN ); check( ); }
};
istream &operator >> ( istream &is, bigInt &tpN ){
    string s;
    is >> s; tpN.clear( );
    for( int i = s.size( ) - 1; i >= 0; --i ) tpN.push_back( s[i] - '0' );
    return is;
}
ostream &operator << ( ostream &os, const bigInt &tpN ){
    if( tpN.empty( ) )    os << 0;
    for( int i = tpN.size( ) - 1; i >= 0; --i )    os << tpN[i];
    return os;
}
bool operator != ( const bigInt &one, const bigInt &another ){
    if( one.size( ) != another.size( ) )    return 1;
    for( int i = one.size( ) - 1; i >= 0; --i ){
        if( one[i] != another[i] )    return 1;
    }
    return 0;
}
bool operator == ( const bigInt &one, const bigInt &another ){
    return ! ( one != another );
}
bool operator < ( const bigInt &one, const bigInt &another ){
    if( one.size( ) != another.size( ) )    return one.size( ) < another.size( );
    for( int i = one.size( ) - 1; i >= 0; --i ){
        if( one[i] != another[i] )    return one[i] < another[i];
    }
    return 0;
}
bool operator > ( const bigInt &one, const bigInt &another ){ return another < one; }
bool operator <= ( const bigInt &one, const bigInt &another ){ return ! (one > another ); }
bool operator >= ( const bigInt &one, const bigInt &another ){ return ! (one < another ); }
bigInt &operator += ( bigInt &one, const bigInt &another ){
    if( one.size( ) < another.size( ) )    one.resize(another.size( ) );
    for( unsigned i = 0; i != another.size( ); ++ i ) one[i] += another[i];
    return one.check( );
}
bigInt operator + ( bigInt one, const bigInt &another ){ return one += another; }
bigInt &operator -= ( bigInt &one, bigInt another ){
    if( one < another )    swap( one, another );
    for( unsigned i = 0; i != another.size( ); one[i] -= another[i], ++ i ){
        if( one[i] < another[i] ){
            unsigned j = i + 1;
            while( ! one[j] ) ++ j;
            while( j > i ){ -- one[j]; one[--j] += 10; }
        }
    }
    return one.check( );
}
bigInt operator - ( bigInt one, const bigInt &another ){ return one -= another; }
bigInt operator * ( const bigInt &one, const bigInt &another ){
    bigInt tpN;
    tpN.assign( one.size( ) + another.size( ) - 1, 0 );
    for( unsigned i = 0; i != one.size( ); ++ i ){
        for( unsigned j = 0; j != another.size( ); ++ j ) tpN[i + j] += one[i] * another[j];
    }
    return tpN.check( );
}
bigInt &operator *= ( bigInt &one, const bigInt &another ){ return one = one * another; }
bigInt divMod( bigInt &one, const bigInt &another ){
    bigInt ans;
    for( int t = one.size( ) - another.size( ); one >= another; -- t ){
        bigInt tpS;
        tpS.assign( t + 1, 0 );
        tpS.back( ) = 1;
        bigInt tpM = another * tpS;
        while( one >= tpM ){ one -= tpM; ans += tpS; }
    }
    return ans;
}
bigInt operator / ( bigInt one, const bigInt &another ){ return divMod(one, another ); }
bigInt &operator /= ( bigInt &one, const bigInt &another ){ return one = one / another; }
bigInt &operator %= ( bigInt &one, const bigInt &another ){ divMod( one, another ); return one; }
bigInt operator % ( bigInt one, const bigInt &another ){ return one %= another; }

struct matrixS{
    int mat[2][2];
    matrixS( int x = 0 ){ memset( mat, x, sizeof( mat ) ); }
    matrixS operator * ( const matrixS &another ) const{
        matrixS res;
        for( int i = 0; i < 2; ++ i ){
            for( int j = 0; j < 2; ++ j ){
                for( int k = 0; k < 2; ++ k )    res.mat[i][j] = add( ( LL )res.mat[i][j], ( LL )mat[i][k] * another.mat[k][j] );
            }
        }
        return res;
    }
} unit, erng;

bigInt N;

void progressBaseInformation( ){
    int unitS[2][2] = { { 1, 0 }, { 0, 1 } };
    memcpy( unit.mat, unitS, sizeof( unitS ) );
    int erngS[2][2] = { { 2, 1 }, { 1, 0 } };
    memcpy( erng.mat, erngS, sizeof( erngS ) );
}

signed main( ){
    ios::sync_with_stdio( 0 ); cin.tie( 0 ); cout.tie( 0 );
    progressBaseInformation( );
    cin >> N; cout << qkpow( erng, unit, N ).mat[1][0] << '\n';
    return 0;
}

不,凉心出题人友好地卡了高精的常数,于是你打开题解,发现 $f_{n}=f_{n\bmod (10^{9}+6)}$,于是你又行了。

$\mathcal{Code}$

#include <cstdio>
#include <cstring>
#include <queue>
#define mod ( 1000000007 )

using namespace std;
typedef long long LL;

template<typename _T>
void read( _T &x ){
    x = 0; char c = getchar( ); _T f = 1;
    while( c < '0' || c > '9' ){ if( c == '-' )    f = -1; c = getchar( ); }
    while( c >= '0' && c <= '9' ){ x = ( ( x << 3 ) + ( x << 1 ) + ( c & 15 ) ) % ( mod - 1 ); c = getchar( ); }
    x *= f;
}

template<typename _T>
void write( _T x ){
    if( x < 0 ){ putchar( '-' ); x = -x; }
    if( x > 9 )    write( x / 10 );
    putchar( x % 10 + '0' );
}

template<typename _T, typename _P>
_T qkpow( _T bas, _T one, _P fur ){
    _T res = one;
    while( fur != 0 ){
        if( fur % 2 == ( _P )1 )    res = bas * res;
        bas = bas * bas;
        fur /= 2;
    }
    return res;
}

template<typename _T>
_T add( _T x, _T y ){ if( y >= mod )    y %= mod; x += y; if( x >= mod )    x -= mod; return x; }

struct matrixS{
    int mat[2][2];
    matrixS( int x = 0 ){ memset( mat, x, sizeof( mat ) ); }
    matrixS operator * ( const matrixS &another ) const{
        matrixS res;
        for( int i = 0; i < 2; ++ i ){
            for( int j = 0; j < 2; ++ j ){
                for( int k = 0; k < 2; ++ k )    res.mat[i][j] = add( ( LL )res.mat[i][j], ( LL )mat[i][k] * another.mat[k][j] );
            }
        }
        return res;
    }
} unit, erng;

LL N;

void progressBaseInformation( ){
    int unitS[2][2] = { { 1, 0 }, { 0, 1 } };
    memcpy( unit.mat, unitS, sizeof( unitS ) );
    int erngS[2][2] = { { 2, 1 }, { 1, 0 } };
    memcpy( erng.mat, erngS, sizeof( erngS ) );
}

signed main( ){
    progressBaseInformation( );
    read( N ); write( qkpow( erng, unit, N ).mat[1][0] ), putchar( '\n' );
    return 0;
}

Description

Link.

这道题 的加强版。

Solution

题解里面大多数都是概率 DP,或者是期望 DP 然后是逆推。甚至不给 DP 的转移式。机房 yyds Reanap 发了一篇逆推的题解,那我就来补一篇正推的期望 DP 的填表法做法。

首先这道题看上去好像可以状压的样子,我们可以设 $f_{i,S}$ 表示当前打了 $i$ 次,敌方的情况是 $S$ 的期望。

不过仔细想一下发现我们只需要知道各种血量的奴隶主有多少即可。

于是我们重新设计 DP 的状态:$f_{s,i,j,k}$ 表示目前打了 $s$ 次,敌方分别有 $i$、$j$、$k$ 个 1hp、2hp、3hp 的奴隶主。

首先我们令 $T=[i+j+k<K]$

那么我们的方程就是:

$$ f_{s,i,j,k}=\begin{cases} f_{s-1,i-1,j,k}+\frac{i}{i+j+k+1},M=1\land i\neq0 \\ f_{s-1,i-1,j,k}+\frac{i}{i+j+k+1},M=2\land i\neq0 \\ f_{s-1,i+1,j-1+T,k}+\frac{j}{i+j+k+1},M=2\land j\neq0 \\ f_{s-1,i-1,j,k}+\frac{i}{i+j+k+1},M=3\land i\neq0 \\ f_{s-1,i+1,j-1,k+T}+\frac{j}{i+j+k+1},M=3\land j\neq0 \\ f_{s-1,i,j+1,k-1+T}+\frac{k}{i+j+k+1},M=3\land k\neq0 \end{cases} $$

这个方程挺好理解的,基本就等于照题意模拟。

然后我们发现转移式中的系数部分和 $f$ 数组没有关系,所以我们可以用矩阵来加速这个东西。

数一数状态数,直接加速直接 T 飞。

有一个矩阵加速常用的 trick,预处理矩阵 2 的幂。

然后取模卡卡常即可。

(代码不保证稳定能过)

#include <cstdio>
#define mod ( 998244353 )

using namespace std;
typedef long long LL;

char buf[1 << 21], *p1 = buf, *p2 = buf;
#define getchar( ) ( p1 == p2 && ( p2 = ( p1 = buf ) + fread( buf, 1, 1 << 21, stdin ), p1 == p2 ) ? EOF : *p1 ++ )

template<typename _T>
void read( _T &x )
{
    x = 0;
    char c = getchar( );
    _T f = 1;
    while( c < '0' || c > '9' )
    {
        if( c == '-' )    f = -1;
        c = getchar( );
    }
    while( c >= '0' && c <= '9' )    x = ( x << 3 ) + ( x << 1 ) + ( c ^ '0' ), c = getchar( );
    x *= f;
}

template<typename _T, typename... Args>
void read( _T &t, Args&... args ) { read( t ), read( args... ); }

template<typename _T>
void write( _T x )
{
    if( x < 0 )    putchar( '-' ), x = -x;
    if( x > 9 )    write( x / 10 );
    putchar( x % 10 + '0' );
}

template<typename _T>
void Add( _T &x, _T y )
{
    if( y >= mod )  y %= mod;
    x += y;
    if( x >= mod )  x -= mod;
}

template<typename _T>
_T square( _T x ) { return x * x; }

const int Maxn = 10 + 5, Maxk = 170 + 5;
int T, M, K, S, Unite[Maxn][Maxn][Maxn];
LL tmp[Maxk], Ans[Maxk], Inv[Maxn];
struct Matrix
{
    LL mat[Maxk][Maxk];

    friend Matrix operator * ( const Matrix &one, const Matrix &another )
    {
        Matrix res;
        for( int i = 0; i <= S + 1; ++ i )
        {
            for( int j = 0; j <= S + 1; ++ j )
            {
                res.mat[i][j] = 0;
                for( int k = 0; k <= S + 1; ++ k )    Add( res.mat[i][j], one.mat[i][k] * another.mat[k][j] );
            }
        }
        return res;
    }
} dp[Maxk];

template<typename _T>
_T qkpow( _T base, _T times )
{
    _T res = 1;
    while( times )
    {
        if( times & 1 )    res = ( LL )res * base % mod;
        base = ( LL )base * base % mod;
        times >>= 1;
    }
    return res;
}

void progressInversions( ) { for( int i = 0; i <= 10; ++ i )    Inv[i] = qkpow( i, mod - 2 ); }
signed main( )
{
    progressInversions( );
    read( T, M, K );
    for( int i = 0; i <= K; ++ i )
    {
        int UpI;
        if( M > 1 )    UpI = K - i;
        else    UpI = 0;
        for( int j = 0; j <= UpI; ++ j )
        {
            int UpJ;
            if( M > 2 )    UpJ = K - i - j;
            else    UpJ = 0;
            for( int k = 0; k <= UpJ; ++ k )    Unite[i][j][k] = ++ S;
        }
    }
    for( int i = 0; i <= K; ++ i )
    {
        int UpI;
        if( M > 1 )    UpI = K - i;
        else    UpI = 0;
        for( int j = 0; j <= UpI; ++ j )
        {
            int UpJ;
            if( M > 2 )    UpJ = K - i - j;
            else    UpJ = 0;
            for( int k = 0; k <= UpJ; ++ k )
            {
                int Add;
                if( i + j + k < K )    Add = 1;
                else    Add = 0;
                if( M == 1 && i )    dp[0].mat[Unite[i][j][k]][Unite[i - 1][j][k]] = ( LL )i * Inv[i + j + k + 1] % mod;
                else if( M == 2 )
                {
                    if( i )    dp[0].mat[Unite[i][j][k]][Unite[i - 1][j][k]] = ( LL )i * Inv[i + j + k + 1] % mod;
                    if( j )    dp[0].mat[Unite[i][j][k]][Unite[i + 1][j - 1 + Add][k]] = ( LL )j * Inv[i + j + k + 1] % mod;
                }
                else if( M == 3 )
                {
                    if( i )    dp[0].mat[Unite[i][j][k]][Unite[i - 1][j][k]] = ( LL )i * Inv[i + j + k + 1] % mod;
                    if( j )    dp[0].mat[Unite[i][j][k]][Unite[i + 1][j - 1][k + Add]] = ( LL )j * Inv[i + j + k + 1] % mod;
                    if( k )    dp[0].mat[Unite[i][j][k]][Unite[i][j + 1][k - 1 + Add]] = ( LL )k * Inv[i + j + k + 1] % mod;
                }
                dp[0].mat[Unite[i][j][k]][Unite[i][j][k]] = dp[0].mat[Unite[i][j][k]][S + 1] = Inv[i + j + k + 1];
            }
        }
    }
    dp[0].mat[S + 1][S + 1] = 1;
    for( int i = 1; i <= 60; ++ i )    dp[i] = square( dp[i - 1] );
    while( T -- > 0 )
    {
        LL N;
        read( N );
        for( int i = 0; i <= S + 1; ++ i )  Ans[i] = 0;
        if( M == 1 )    Ans[Unite[1][0][0]] = 1;
        else if( M == 2 )    Ans[Unite[0][1][0]] = 1;
        else    Ans[Unite[0][0][1]] = 1;
        for( int i = 0; i <= 60; ++ i )
        {
            if( ( N >> i ) & 1 )
            {
                for( int j = 0; j <= S + 1; ++ j )
                {
                    tmp[j] = 0;
                    for( int k = 0; k <= S + 1; ++ k )    Add( tmp[j], Ans[k] * dp[i].mat[k][j] );
                }
                for( int j = 0; j <= S + 1; ++ j )    Ans[j] = tmp[j];
            }
        }
        write( Ans[S + 1] ), putchar( '\n' );
    }
    return 0;
}