板板。
#include<bits/stdc++.h> #define ll long long #define N 22 #define P 998244353 ll t,p,q,r,l; struct Po{ ll cntu,cntr,sumi,sums,sqrs,prod; Po(){cntu = cntr = sumi = sums = sqrs = prod = 0;} Po operator + (Po b) { Po c; c.cntu=(cntu+b.cntu)%P,c.cntr=(cntr+b.cntr)%P; c.sumi=(sumi+b.sumi+cntr*b.cntr)%P; c.sums=(sums+b.sums+cntu*b.cntr)%P; c.sqrs=(sqrs+b.sqrs+((cntu*cntu)%P)*b.cntr+(2*cntu*b.sums)%P)%P; c.prod=((prod+b.prod+((cntu*cntr)%P)*b.cntr)%P+cntu*b.sumi+cntr*b.sums)%P; return c; } }nu,nr,ans; inline Po pow(Po a,ll k){ Po res; while(k){ if(k & 1){res = res + a;} a = a + a; k >>= 1; } return res; } inline ll div(ll a,ll b,ll c,ll d){return ((long double)1.0 * a * b + c) / d;} inline Po solve(ll p,ll q,ll r,ll l,Po a,Po b){ if(!l)return Po(); if(p >= q)return solve(p % q,q,r,l,a,pow(a,p / q) + b); ll m = div(l,p,r,q); if(!m)return pow(b,l); ll cnt = l - div(q,m,-r - 1,p); return pow(b,(q - r - 1) / p) + a + solve(q,p,(q - r - 1) % p,m - 1,b,a) + pow(b,cnt); } int main(){ scanf("%lld",&t); while(t -- ){ scanf("%lld%lld%lld%lld",&l,&p,&r,&q); nu.cntu = 1,nu.cntr = nu.sumi = nu.sums = nu.sqrs = nu.prod = 0; nr.cntu = nr.sums = nr.sqrs = nr.prod = 0,nr.sumi = nr.cntr = 1; ans = pow(nu,r / q) + solve(p,q,r % q,l,nu,nr); printf("%lld %lld %lld\n",(ans.sums+r/q)%P,(ans.sqrs+((r/q)%P)*((r/q)%P))%P,ans.prod); } }
把其看做\(\sum a^xb^y\)类型。
则可写作矩阵:
\(\begin{bmatrix}\sum_x a^xb^y\\a^xb^y\end{bmatrix}\)\(\begin{bmatrix}1&0\\0&b\end{bmatrix}\)\(\begin{bmatrix}1&a\\0&1\end{bmatrix}\)
那么把\(a,b\)写作矩阵形式也可以
先判断一个数是否是质数,使用Miller Rabin测试,否则用Pollard Rho算法 找到一个因数,递归操作\(n / p\)和\(p\)。
#include <bits/stdc++.h> using namespace std; typedef long long ll; int t; long long max_factor, n; long long gcd(long long a, long long b) { if (b == 0) return a; return gcd(b, a % b); } long long quick_pow(long long x, long long p, long long mod) { //快速幂 long long ans = 1; while (p) { if (p & 1) ans = (__int128)ans * x % mod; x = (__int128)x * x % mod; p >>= 1; } return ans; } bool Miller_Rabin(long long p) { //判断素数 if (p < 2) return 0; if (p == 2) return 1; if (p == 3) return 1; long long d = p - 1, r = 0; while (!(d & 1)) ++r, d >>= 1; //将d处理为奇数 for (long long k = 0; k < 10; ++k) { long long a = rand() % (p - 2) + 2; long long x = quick_pow(a, d, p); if (x == 1 || x == p - 1) continue; for (int i = 0; i < r - 1; ++i) { x = (__int128)x * x % p; if (x == p - 1) break; } if (x != p - 1) return 0; } return 1; } long long Pollard_Rho(long long x) { long long s = 0, t = 0; long long c = (long long)rand() % (x - 1) + 1; int step = 0, goal = 1; long long val = 1; for (goal = 1;; goal *= 2, s = t, val = 1) { //倍增优化 for (step = 1; step <= goal; ++step) { t = ((__int128)t * t + c) % x; val = (__int128)val * abs(t - s) % x; if ((step % 127) == 0) { long long d = gcd(val, x); if (d > 1) return d; } } long long d = gcd(val, x); if (d > 1) return d; } } void fac(long long x) { if (x <= max_factor || x < 2) return; if (Miller_Rabin(x)) { //如果x为质数 max_factor = max(max_factor, x); //更新答案 return; } long long p = x; while (p >= x) p = Pollard_Rho(x); //使用该算法 while ((x % p) == 0) x /= p; fac(x), fac(p); //继续向下分解x和p } int main() { scanf("%d", &t); while (t--) { srand((unsigned)time(NULL)); max_factor = 0; scanf("%lld", &n); fac(n); if (max_factor == n) //最大的质因数即自己 printf("Prime\n"); else printf("%lld\n", max_factor); } return 0; }