[noi simulation] pendulum (linear algebra, Du Jiao sieve)

Topic

6s , 1024mb

I am a XYX, I'm good at swinging .

I saw one when I was putting it on $n\times n$ Matrix $A$ ：
$$A_{i,j}=\{\begin{array}{ll}1& i=j\\ 0& i\ne j\wedge i|j\\ C& \mathrm{o}\mathrm{t}\mathrm{h}\mathrm{e}\mathrm{r}\mathrm{w}\mathrm{i}\mathrm{s}\mathrm{e}\end{array}$$

Now I want to know $A$ The determinant of , Because I put , So the answer just needs to be right $998244353$ modulus .

$1\le n\le 10^{11},0\le C<998244353$ .

Examples ：

6 2$\to$998244350

99 4$\to$714389845

10000000000 10$\to$82177551

Answer key

This problem tells us , It is not necessary to construct upper triangular or lower triangular matrix to calculate determinant .

We can also construct Heisenberg matrix relatively easily .

This is the matrix $A$ ：
$$\left[\begin{array}{cccccc}1& 0& 0& 0& 0& \dots \\ C& 1& C& 0& C& \dots \\ C& C& 1& C& C& \dots \\ C& C& C& 1& C& \dots \\ C& C& C& C& 1& \dots \\ \dots & \dots & \dots & \dots & \dots & \ddots \end{array}\right]$$

Let's put the matrix $A$ Subtract the upper line from each line （ Start from below , Subtract the previous line from each line , The first line doesn't move ）, We can get a Shanghai senberg matrix ：
$$\left[\begin{array}{cccccc}1& 0& 0& 0& 0& \dots \\ C-1& 1& C& 0& C& \dots \\ 0& C-1& 1-C& C& 0& \dots \\ 0& 0& C-1& 1-C& 0& \dots \\ 0& 0& 0& C-1& 1-C& \dots \\ \dots & \dots & \dots & \dots & \dots & \ddots \end{array}\right]$$

Each permutation ring of the matrix is a counter clockwise ring with continuous labels , Be similar to $i\to j\to j-1\to ...\to i+1$ Of . Differential obtainable , Each multiple relationship in the upper right $i|j$ , Will give $A_{i,j}$ contribution $-C$ , to $A_{i+1,j}$ contribution $C$. We multiply each position of the matrix by $\frac{1}{1-C}$ , Then make $f_i$ Before retention $i$ That's ok $i$ Determinant of column , that
$$f_i=f_{i-1}+\frac{C}{1-C}\sum _{j|i,j<i}(f_j-f_{j-1})$$

set up $g_i=f_i-f_{i-1}$ , that $g_i=\frac{C}{1-C}{\sum }_{j|i,j<i}{g}_j$ .

So what we want is $g_i$ And $f_n$ , The data range should be a sub linear sieve .

Make $h_1=-1,h_i=\frac{C}{1-C}$ , that （ Dirichlet convolution ） $h\ast g=-e$ , So you can Dujiao sieve （ Duchenne sieve is not restricted to be a product function ）.
$$\begin{gathered} \sum _{i=1}^n\sum _{j|i}{h}_j{g}_{i/j}=\sum _j^n{h}_j\sum _i^{n/j}{g}_i=-1 \\ \Rightarrow f_n=1+\sum _{j=2}^n{h}_j{f}_{n/j} \end{gathered}$$

We preprocess $n^{2/3}$ Inside $f_i$ , So the complexity of Du Jiao sieve is $O(n^{2/3})$ .

But preprocessing can easily bring a $\log$ , So we need to optimize .

Make $x=\prod p_i^{w_i}$ , We put $p_i$ Change from small to large $2,3,5,7,...$ obtain $y=\prod p{{}^{\prime }}_i^{w_i}$ , Easy to find $g_x=g_y$ , because $g$ Only follow $w_i$ On the resettable set of . We use Euler sieve to find the maximum prime factor and the number of prime factors for each number , And then recursively calculate each $x$ Corresponding $y$ . Violence finds something fundamentally different $y$ Only a thousand at most , We can each $O(n^{1/3})$ Find out $g$ .

Be careful , Because every position of the matrix is multiplied by $\frac{1}{1-C}$ , And because the upper left corner of the matrix is not standard , So the final answer is
$$f_n\cdot (1-C{)}^{n-1}$$

The sum $O(n^{2/3})$ .

CODE

#include<map>
#include<set>
#include<cmath>
#include<ctime>
#include<queue>
#include<stack>
#include<random>
#include<bitset>
#include<vector>
#include<cstdio>
#include<cstring>
#include<iostream>
#include<algorithm>
#include<unordered_map>
#pragma GCC optimize(2)
using namespace std;
#define MAXN 10000005
#define LL long long
#define ULL unsigned long long
#define ENDL putchar('\n')
#define DB double
#define lowbit(x) (-(x) & (x))
#define FI first
#define SE second
#define PR pair<int,int>
int xchar() {
    
	static const int maxn = 1000000;
	static char b[maxn];
	static int pos = 0,len = 0;
	if(pos == len) pos = 0,len = fread(b,1,maxn,stdin);
	if(pos == len) return -1;
	return b[pos ++];
}
// #define getchar() xchar()
LL read() {
    
	LL f = 1,x = 0;int s = getchar();
	while(s < '0' || s > '9') {
    if(s<0)return -1;if(s=='-')f=-f;s = getchar();}
	while(s >= '0' && s <= '9') {
    x = (x<<1) + (x<<3) + (s^48);s = getchar();}
	return f*x;
}
void putpos(LL x) {
    if(!x)return ;putpos(x/10);putchar((x%10)^48);}
void putnum(LL x) {
    
	if(!x) {
    putchar('0');return ;}
	if(x<0) putchar('-'),x = -x;
	return putpos(x);
}
void AIput(LL x,int c) {
    putnum(x);putchar(c);}

const int MOD = 998244353;
int n,m,s,o,k;
LL N;
inline int qkpow(int a,int b) {
    
	int res = 1;
	while(b > 0) {
    
		if(b & 1) res = res *1ll* a % MOD;
		a = a *1ll* a % MOD; b >>= 1;
	} return res;
}
int V;
int sg[MAXN],lw[MAXN],ct[MAXN],hb[MAXN];
int p[MAXN],cnt; bool f[MAXN];
void sieve(int n) {
    
	sg[1] = 1; lw[1] = 1; int nm = 0;
	for(int i = 2;i <= n;i ++) {
    
		if(!f[i]) p[++ cnt] = i,lw[i] = 2,hb[i] = i,ct[i] = 1;
		if(lw[i] == i) {
    
			for(int j = 1;j * j <= i;j ++) {
    
				if(i % j == 0) {
    
					(sg[i] += sg[j]) %= MOD;
					if(i/j > j && j > 1) {
    
						(sg[i] += sg[i/j]) %= MOD;
					}
				}
			}
			sg[i] = sg[i] *1ll* V % MOD;
		}
		else sg[i] = sg[lw[i]];
		for(int j = 1;j <= cnt && (nm=i*p[j]) <= n;j ++) {
    
			f[nm] = 1; hb[nm] = hb[i];
			if(i % p[j] == 0) {
    
				ct[nm] = ct[i];
				lw[nm] = lw[nm/hb[nm]] * p[ct[nm]];
				break;
			}
			ct[nm] = ct[i] + 1;
			lw[nm] = lw[nm/hb[nm]] * p[ct[nm]];
		}
	}
	for(int i = 2;i <= n;i ++) {
    
		(sg[i] += sg[i-1]) %= MOD;
	}
	return ;
}
int sh(LL x) {
    
	if(x<1) return 0;
	return ((x-1)%MOD*V%MOD +MOD-1) % MOD;
}
int s1[1000005];
int Sg(LL x) {
    
	if(x <= 10000000) return sg[x];
	int &rs = s1[N/x];
	if(rs >= 0) return rs;
	rs = 1;
	for(LL i = 2;i <= x;i ++) {
    
		LL r = x/(x/i);
		rs += (r-i+1)%MOD*V%MOD *1ll* Sg(x/i) % MOD;
		if(rs >= MOD) rs -= MOD;
		i = r;
	} return rs;
}
int main() {
    
	freopen("bigben.in","r",stdin);
	freopen("bigben.out","w",stdout);
	N = read(); int C = read();
	if(C <= 1) return AIput(N<=2,'\n'),0;
	V = C *1ll* qkpow(MOD+1-C,MOD-2) % MOD;
	sieve(10000000);
	memset(s1,-1,sizeof(s1));
	int as = Sg(N) *1ll* qkpow(MOD+1-C,(N-1)%(MOD-1)) % MOD;
	AIput(as,'\n');
	return 0;
}