[SDOI2013] 随机数生成器

description

solution

肯定是非常想找一个通项公式来表示第$n$个数的

依据形式，考虑化成等比数列

$x_{i+1}+k=a(x_i+k)=a\cdot x_i+b+t\Rightarrow k=\frac{b}{a-1}$

$\Rightarrow x_i+\frac{b}{a-1}=a^{i-1}(x_1+\frac{b}{a-1})(\mathrm{m}\mathrm{o}\mathrm{d}p)\Rightarrow a^{i-1}\equiv (x_i+b\cdot (a-1{)}^{-1})\cdot (x_1+b\cdot (a-1{)}^{-1}{)}^{-1}(\mathrm{m}\mathrm{o}\mathrm{d}p)$

$x_i$就是题目中的$t$，所以式子只有$a^{i-1}$这个指数未知，bsgs求解最小指数（$p$是质数）

但要先特判$a=0/1$两种情况，显然此情况下等比数列式子是不成立的

• a=0

  $x_i=b$

• a=1

  $x_i\equiv x_1+b(i-1)(\mathrm{m}\mathrm{o}\mathrm{d}p)\iff (t-x_1)b^{-1}\equiv i-1(\mathrm{m}\mathrm{o}\mathrm{d}p)$

code

#include <map>
#include <cmath>
#include <cstdio>
using namespace std;
#define int long long
map < int, int > mp;int bsgs( int a, int r, int p ) {if( a % p == 0 ) return -1;mp.clear();int m = ceil( sqrt( p * 1.0 ) );int IS = 1;for( int i = 1;i <= m;i ++ ) {IS = IS * a % p;mp[IS * r % p] = i;}int SI = 1;for( int i = 1;i <= m;i ++ ) {SI = SI * IS % p;if( mp.count( SI ) ) return i * m - mp[SI] + 1;}return -1;
}int T, p, a, b, x1, t;int qkpow( int x, int y ) {int ans = 1;while( y ) {if( y & 1 ) ans = ans * x % p;x = x * x % p;y >>= 1;}return ans;
}int inv( int x ) {return qkpow( x, p - 2 );
}signed main() {scanf( "%lld", &T );while( T -- ) {scanf( "%lld %lld %lld %lld %lld", &p, &a, &b, &x1, &t );if( x1 == t ) {printf( "1\n" );continue;}if( a == 0 ) {if( x1 == t ) printf( "1\n" );else if( b == t ) printf( "2\n" );else printf( "-1\n" );continue;}if( a == 1 ) {if( b == 0 ) printf( "-1\n" );else printf( "%lld\n", ( ( t - x1 ) % p + p ) % p * inv( b ) % p + 1 );continue;}printf( "%lld\n", bsgs( a, ( t + b * inv( a - 1 ) % p ) * inv( x1 + b * inv( a - 1 ) % p ) % p, p ) );}	return 0;
}

[BZOJ4128]Matrix

description

solution

因为保证了$A$矩阵有逆，所以可以选择求矩阵的逆

但其实还是选择两边同乘矩阵要好一些

BSGS模板，有优质代码是写了非常迅速的hash

但是偶直接map里面塞矩阵，重载一下map的小于等于比较挺不戳的！

code

#include <map>
#include <cmath>
#include <cstdio>
#include <cstring>
using namespace std;
int n, mod;
struct matrix {int c[71][71];void clear() {memset( c, 0, sizeof( c ) );}matrix() {clear();}void init() {for( int i = 1;i <= n;i ++ ) {for( int j = 1;j <= n;j ++ )c[i][j] = 0;c[i][i] = 1;}}matrix operator * ( const matrix &t ) const {matrix ans;for( int i = 1;i <= n;i ++ )for( int j = 1;j <= n;j ++ )for( int k = 1;k <= n;k ++ )ans.c[i][j] = ( ans.c[i][j] + c[i][k] * t.c[k][j] ) % mod;return ans;}friend bool operator == ( const matrix &s, const matrix &t ) {for( int i = 1;i <= n;i ++ )for( int j = 1;j <= n;j ++ )if( s.c[i][j] ^ t.c[i][j] ) return 0;return 1;}friend bool operator < ( const matrix &s, const matrix &t ) {for( int i = 1;i <= n;i ++ )for( int j = 1;j <= n;j ++ )if( s.c[i][j] ^ t.c[i][j] ) return s.c[i][j] < t.c[i][j];return 0;}
}A, B, g, h;
map < matrix, int > mp;int main() {scanf( "%d %d", &n, &mod );for( int i = 1;i <= n;i ++ )for( int j = 1;j <= n;j ++ )scanf( "%d", &A.c[i][j] );for( int i = 1;i <= n;i ++ )for( int j = 1;j <= n;j ++ )scanf( "%d", &B.c[i][j] );int m = ceil( sqrt( mod * 1.0 ) );g.init();for( int i = 1;i <= m;i ++ ) {g = g * A;h = g * B;mp[h] = i;}h.init();for( int i = 1;i <= m;i ++ ) {h = h * g;if( mp.count( h ) ) return ! printf( "%d\n", i * m - mp[h] );}return 0;
}

CF1106F Lunar New Year and a Recursive Sequence

problem

solution

$f_i={\prod }_{j=1}^k{f}_{i-j}^{b_j}(\mathrm{m}\mathrm{o}\mathrm{d}p)$

$p$是非常特殊的质数，其原根为$3$，可以用原根改写上式，令$3^{g_i}=f_i$

$g_i={\sum }_{j=1}^k{g}_{i-j}\cdot b_j(\mathrm{m}\mathrm{o}\mathrm{d}\phi (p))$

此时的$g$可以用矩阵加速，初始化$f_i=0(i<k),f_k=1$，转移$n-k$次变成$[f_{n-k+1},...,f_n]$
$$\left[\begin{array}{c}\underset{k}{\underbrace{h_{i-k+1},...,h_i}}\end{array}\right]\times \left[\begin{array}{ccccc}0& 0& ...& 0& b_k\\ 1& 0& ...& 0& b_{k-1}\\ 0& 1& ...& 0& b_{k-2}\\ ...& ...& ...& ...& ...\\ 0& 0& ...& 0& b_2\\ 0& 0& ...& 1& b_1\end{array}\right]=\left[\begin{array}{c}\underset{k}{\underbrace{h_{i+k+2},...,h_{i+1}}}\end{array}\right]$$
同时，有$f_k^t\equiv f_n(\mathrm{m}\mathrm{o}\mathrm{d}p),3^{g_k}=f_k$，这里的$t$矩等于加速矩阵自乘后的$h_k$

$f_k^T=3^{g_k\cdot T}\equiv f_n=3^{g_n}(\mathrm{m}\mathrm{o}\mathrm{d}p)$

对于$3^{g_n}\equiv f_n(\mathrm{m}\mathrm{o}\mathrm{d}p)$的$g_n$未知，用BSGS求得

对于指数用exgcd计算，$g_k\times t\equiv g_n(\mathrm{m}\mathrm{o}\mathrm{d}p-1)\iff t\times g_k-(p-1)\times y=g_n$

未知数为$g_k,y$，求出后$f_k=3^{g_k}$

code

#include <map>
#include <cmath>
#include <cstdio>
#include <cstring>
using namespace std;
#define int long long
#define mod 998244353
struct matrix {int n, m;int c[101][101];matrix() {memset( c, 0, sizeof( c ) );}matrix operator * ( matrix &t ) const {matrix ans;ans.n = n, ans.m = t.m;for( int i = 1;i <= n;i ++ )for( int j = 1;j <= t.m;j ++ )for( int k = 1;k <= m;k ++ )ans.c[i][j] = ( ans.c[i][j] + c[i][k] * t.c[k][j] ) % ( mod - 1 );return ans;}
}h, t;matrix qkpow( matrix x, int y ) {matrix ans;ans.n = ans.m = x.n;for( int i = 1;i <= ans.n;i ++ )ans.c[i][i] = 1;while( y ) {if( y & 1 ) ans = ans * x;x = x * x;y >>= 1;}return ans;
}map < int, int > mp;int bsgs( int a, int r ) {int m = ceil( sqrt( mod * 1.0 ) );int IS = 1;for( int i = 1;i <= m;i ++ ) {IS = IS * a % mod;mp[IS * r % mod] = i;}int MS = 1;for( int i = 1;i <= m;i ++ ) {MS = MS * IS % mod;if( mp.count( MS ) )return i * m - mp[MS];}return -1;
}int exgcd( int a, int b, int &x, int &y ) {if( ! b ) {x = 1, y = 0;return a;}else {int d = exgcd( b, a % b, y, x );y -= a / b * x;return d;}
}int qkpow( int x, int y ) {int ans = 1;while( y ) {if( y & 1 ) ans = ans * x % mod;x = x * x % mod;y >>= 1;}return ans;
}signed main() {int k, n, fn;scanf( "%lld", &k );h.n = h.m = t.m = k, t.n = t.c[1][k] = 1;for( int i = 2;i <= k;i ++ )h.c[i][i - 1] = 1;for( int i = k;i;i -- )scanf( "%lld", &h.c[i][k] );scanf( "%lld %lld", &n, &fn );h = qkpow( h, n - k );t = t * h;int g = bsgs( 3, fn );int x, y, d = exgcd( t.c[1][k], mod - 1, x, y );if( g % d ) return ! printf( "-1\n" );else {x = ( x * ( g / d ) % ( mod - 1 ) + ( mod - 1 ) ) % ( mod - 1 );printf( "%lld\n", qkpow( 3, x ) );}return 0;
}

Fermat’s Last Theorem

description

solution

找出$p$的原根$g$，令$x=g^a,y=g^b,z=g^c$
则$g^{an}+g^{bn}\equiv g^{cn}(\mathrm{m}\mathrm{o}\mathrm{d}p)$
令$r=g^n\Rightarrow r^a+r^b\equiv r^c(\mathrm{m}\mathrm{o}\mathrm{d}p)$
假设$a\le b$

• $a\le c\Rightarrow 1+r^{b-a}\equiv r^{c-a}(\mathrm{m}\mathrm{o}\mathrm{d}p)$
• $a>c\Rightarrow r^{a-c}+r^{b-a}\equiv 1(\mathrm{m}\mathrm{o}\mathrm{d}p)$

枚举$r^k$，出现循环的时候就break
出现$1+r^{k_1}\equiv r^{k_2}/r^{k_1}+r^{k_2}\equiv 1(\mathrm{m}\mathrm{o}\mathrm{d}p)$就是有解了

code

#include <cstdio>
#define int long long
#define maxn 1000005
int d[maxn], vis[maxn], mi[maxn];int qkpow( int x, int y, int mod ) {int ans = 1;while( y ) {if( y & 1 ) ans = ans * x % mod;x = x * x % mod;y >>= 1;}return ans;
}int root( int p ) {int cnt = 0, phi = p - 1;for( int i = 2;i * i <= phi;i ++ )if( phi % i == 0 ) {d[++ cnt] = i;if( i * i != phi ) d[++ cnt] = phi / i;}for( int g = 1, i;;g ++ ) {for( i = 1;i <= cnt;i ++ )if( qkpow( g, d[i], p ) == 1 )break;if( i > cnt ) return g;}
}signed main() {int T, n, p;scanf( "%lld", &T );for( int Case = 1;Case <= T;Case ++ ) {scanf( "%lld %lld", &n, &p );int g = root( p ), r = qkpow( g, n, p );int t = 1, x = -1, y, z;for( int i = 0;i < p;i ++ ) {if( i ) t = t * r % p;if( vis[t] == Case ) break;vis[t] = Case, mi[t] = i;int MS = 1 + t;if( MS >= p ) MS = 0;int IS = 1 - t + p;if( IS >= p ) IS = 0;if( vis[MS] == Case ) {x = 1, y = qkpow( g, i, p ), z = qkpow( g, mi[MS], p );break;}if( vis[IS] == Case ) {x = qkpow( g, i, p ), y = qkpow( g, mi[IS], p ), z = 1;break;}}if( ~ x ) printf( "%lld %lld %lld\n", x, y, z );else printf( "-1\n" );}return 0;
}