注意到 $[1,m]$ 都不与 $m!$ 互质，其他不与其互质的数都是 $[1,m]$ 的倍数，考虑到只数出 $[1,n!]$ 有多少个 $[1,m]$ 以内质数的倍数。这玩意和欧拉函数定义式的推导基本一直，读者自行蓝书。

#include<bits/stdc++.h>
using namespace std;
#define int long long
const int N=1e7+12;
int T,mod,cnt;
int mul[N],pri[N/10],ny[N/10];
struct node{
	int n,m,id;
}q[11010];int ans[11010];
bitset<N> vis;
int qpow(int a,int b){
	int res=1;
	while(b){
		if(b&1) res=(res*a)%mod;
		b>>=1;a=(a*a)%mod;
	}return res;
}
bool cmp(node A,node B){
	return A.m<B.m;
}
signed main(){
	ios::sync_with_stdio(0);
	cin.tie(0), cout.tie(0);
	cin >> T >> mod; 
	mul[0]=1; 
	for(int i=1;i<=N-10;++i){ //预处理阶乘 
		mul[i]=(mul[i-1]*i)%mod;
	}
	for(int i=2;i<=N-10;++i){//筛质数 
		if(vis[i]==0){
			pri[++cnt]=i;
		}
		for(int j=1;j<=cnt&&pri[j]*i<=N-10;++j){
			vis[pri[j]*i]=1;
			if(i%pri[j]==0) break;
		}
	}
	for(int i=1;i<=cnt;++i){//对质数处理逆元 
		ny[i]=qpow(pri[i],mod-2);
	}
	for(int i=1;i<=T;++i){//离线处理 
		cin >> q[i].n >> q[i].m;
		q[i].id=i;
	}sort(q+1,q+T+1,cmp);
	int j=1,tmp=1;
	for(int i=1;i<=cnt&&j<=T;++i){
		while(pri[i]>q[j].m){
//			cout << q[j].id <<' '<< q[j].n <<' '<<q[j].m<<" "<<'\n';
			ans[q[j].id]=(mul[q[j].n]*tmp)%mod;
			j++;
			if(j>T) break;
		}
		tmp=(tmp)*((1-ny[i]+mod)%mod)%mod;
	}
	for(int i=1;i<=T;i++){
		cout << ans[i] << '\n';
	}
	return 0;
}
/*

*/

注意到，$[1,m!]$跟$m!$的互质关系与$[km!+1,km!+m!]$的相同 因此可以将$n!$拆成$m!*(n!/m!)$

答案就是$(n!/m!)*phi(m!)$ 把询问离线然后随便做就行

#include<bits/stdc++.h>
//#define int long long 
#define ll long long
//#define big __int128
#define pii pair<int,int>
#define F first
#define S second
#define mkp make_pair 
using namespace std;
const ll inf=1e18;
const ll N=1e7,E=1e4;
bitset<N+E>pvis;
int pri[N/10+E],ptot,mbs[N+E];
ll mod;
struct Qry{
	int n,m,id;
	ll ans;
}qry[10100];
bool cmpn(Qry x,Qry y){
	return x.n<y.n;
}bool cmpm(Qry x,Qry y){
	return x.m<y.m;
}bool cmpid(Qry x,Qry y){
	return x.id<y.id;
}
ll qpow(ll x,ll y){
	ll sum=1;
	while(y){
		if(y&1) sum=sum*x%mod;
		x=x*x%mod;
		y>>=1;
	}return sum;
}
ll J[10100000],Inv[10100000],T;
signed main() {
	ios::sync_with_stdio(0);
	cin.tie(0);
//	freopen("data.in","r",stdin);
//	freopen("my.out","w",stdout);
//	system("fc my.out ex.ans");return 0;
	cin>>T>>mod;
	J[0]=1;int M=min(N,mod-1);
	for(int i=1;i<=M;i++) J[i]=J[i-1]*i%mod;
	Inv[M]=qpow(J[M],mod-2);
	for(int i=M;i>=1;i--)Inv[i-1]=Inv[i]*i%mod;
	for(int i=1;i<=T;i++){
		cin>>qry[i].n>>qry[i].m;qry[i].id=i;
	}
	sort(qry+1,qry+T+1,cmpm);
	int p=1;
	ll phi=1;
	pvis[1]=1;
	for(int i=1;i<=N;i++){
		if(pvis[i]==0){
			phi=phi*(i-1)%mod;
			pri[++ptot]=i;
		}else phi=phi*i%mod;
		for(int j=1;j<=ptot&&pri[j]*i<=N;j++){
			pvis[i*pri[j]]=1;
			if(i%pri[j]==0) break;
		}
		while(qry[p].m==i){
			int qn=qry[p].n,qm=qry[p].m;
			if((qm+mod-1)/mod*mod<=qn) qry[p].ans=0;
			else qry[p].ans=J[qn%mod]*Inv[qm%mod]%mod*phi%mod;
			p++;
		}
	}
	sort(qry+1,qry+T+1,cmpid);
	for(int i=1;i<=T;i++){
		cout<<qry[i].ans<<"\n";
	}
	return 0;
}

要找[1,n!]中与m!互质的数的个数，而与m!互质必须满足不包含任意属于[1,m]的质数（因为m!是由[1,m]中所有质数乘来乘去组成的）

那么（这里反过来想，也可以认为是贡献法）对于属于[1,m]的质数p，[1,n!]中与p不互质的有n!/p个

那么对于一堆p，最后互质的就有 $n!-n!/p_{1}-n!/p_{2}......+容斥$ （因为 $p_{i} 、p_{j}$ 可能会重复，所以还要容斥一下）

这时发现上式与欧拉函数定义式相似，所以容斥方法也一样，-> $n!*(1-1/p_{1})*(1-1/p_{2}).....$

然后就好了，但是单独处理O(nT)，要优化一下，这里显然选择离线，可以对于n或m排个序，预处理n!或p[1]* p[2]* ...* p[m]，可以做到O(n+m/logm*logP) 复杂度O(n)

注意，当n>=P&&m>=P时，n!中 乘的P 会和 <=m的质数中 除的P 相互消掉，这时如果求逆元会使答案为0，而答案可能不为0 所以要特殊处理一下

详见此样例：

1 3
4 3

ans=2

#include<bits/stdc++.h>
using namespace std;
int mb,ct,prime[10000007];
bitset<10000007> not_prime;
void euler(){
	for(int i=2;i<=mb;i++){
		if(!not_prime[i]){
			prime[++ct]=i;
		}
		for(int j=1;j<=ct&&prime[j]*i<=mb;j++){
			not_prime[i*prime[j]]=1;
			if(i%prime[j]==0) break;
		}
	}
}
int P;
struct node{
	int id,a,b;
}s[10004];
bool cmp(node x,node y){
	return x.b<y.b;
}
long long fpw(long long a,int b){
	long long ans=1;
	while(b){
		if(b&1) ans=ans*a%P;
		a=a*a%P;b>>=1;
	}return ans;
}
int T,as[10004],ma;
long long ans1=1,jc[10000007];
int main()
{
	scanf("%d%d",&T,&P);
	for(int i=1;i<=T;i++){
		scanf("%d%d",&s[i].a,&s[i].b);s[i].id=i;
		ma=max(ma,s[i].a);mb=max(s[i].b,mb);
	}
	sort(s+1,s+T+1,cmp);
	euler();
	jc[1]=1;
    for(int i=2;i<=min(ma,2*P);i++){//模P下，P需特殊处理，但2P不是质数不会影响，>=2P就是0了
		jc[i]=jc[i-1];
		if(i!=P) jc[i]=jc[i]*i%P;
	}
	for(int t=1,j=1;t<=T;t++){
		int n=s[t].a,m=s[t].b;
		while(j<=ct&&prime[j]<=m){
			ans1=ans1*(prime[j]-1)%P;
			if(prime[j]!=P){//特判P
				ans1=ans1*fpw(prime[j],P-2)%P;
			}j++;
		}
		if(n>=P&&m>=P||n<P&&m<P){//P会消掉||没有P
			as[s[t].id]=ans1*jc[n]%P;
		}
		else as[s[t].id]=0;//n或m中包含P
	}
	for(int i=1;i<=T;i++){
		printf("%d\n",as[i]);
	}
	return 0;
}