HDU 6169 Senior PanⅡ（数论+dp）

Description

给出一个区间$[L,R]$，问该区间中所有以$K$作为最小因子（大于$1$的）的数字之和

Input

第一行输入一整数$T$表示用例组数，每组用例输入三个整数$L,R,K$$(1\le L\le R\le 10^{11},2\le K\le 10^{11})$

Output

对于每组用例，输出答案，结果模$10^9+7$

Sample Input

2
1 20 5
2 6 3

Sample Output

Case #1: 5
Case #2: 3

Solution

显然$k$需要是素数，否则不存在数字以$K$为最小因子，且如果$K$是$N$的因子，那么$N/K$也是$N$的因子，进而有$K\le N/K$，故$K\le \sqrt{N}$，即如果$K>\sqrt{10^{11}}\approx 320000$时，至多有一个数字（即$K$本身）以$K$为最小因子，当$K\le 320000$时，只需解决$1$~$N$中以$K$为最小素因子的数字之和即可

考虑$dp$，以$dp[i][j]$表示前$i$个数字中没有以前$j$个素数为最小素因子的数字之和，假设$K$是第$pos+1$个因子，如果一个数字$X$以$K$为最小素因子，那么$X/K$必然不以前$pos$个素数为最小素因子，反之，如果$Y$不以前$pos$个素因子为最小素因子，则$YK$必然以$K$为最小素因子，那么我们只要求$1$~$\lfloor\frac{N}{K}\rfloor$中没有以前$pos$个素数为最小素因子的数字之和，乘上$K$即为$1$~$N$中以$K$为最小素因子的数字之和，且同理可以得到转移方程，前$i$个数字中不以前$j$个素数为最小素因子的数字之和$=$前$i$个数字中不以前$j-1$个素数为最小素因子的数字之和$-$前$i$个数字中以第$j$个素数为最小素因子的数字之和，即$dp[i][j]=dp[i][j-1]-p_j\cdot dp[\lfloor\frac{i}{p_j}\rfloor][j-1]$

Code

#include<cstdio>
#include<iostream>
#include<cstring>
#include<algorithm>
#include<cmath>
#include<vector>
#include<queue>
#include<map>
#include<set>
#include<ctime>
using namespace std;
typedef long long ll;
#define M 320000
#define mod 1000000007
#define inv2 500000004
#define A 8000
#define B 1000
int res=0,p[M+5],mark[M+5],dp[A+5][B+5];
int ins(int x,int y)
{
    return x+y>=mod?x+y-mod:x+y;
}
int des(int x,int y)
{
    return x-y<0?x-y+mod:x-y;
}
void init()
{
    for(int i=2;i<=M;i++)
        if(!mark[i])
        {
            p[++res]=i;
            for(int j=2*i;j<=M;j+=i)mark[j]=1;
        }
    dp[0][0]=0;
    for(int i=1;i<=A;i++)
    {
        dp[i][0]=ins(dp[i-1][0],i);
        for(int j=1;j<=B;j++)
            dp[i][j]=des(dp[i][j-1],(ll)p[j]*dp[i/p[j]][j-1]%mod);
    }
}
bool check(ll n)
{
    if(n<=M)return mark[n]^1;
    for(int i=1;i<=res&&(ll)p[i]*p[i]<=n;i++)
        if(n%p[i]==0)return 0;
    return 1;
}
int Solve(ll n,int m)
{
    if(n<2)return n;
    if(m==0)
    {
        n%=mod;
        return n*(n+1)%mod*inv2%mod;
    }
    if(n<=A&&m<=B)return dp[n][m];
    if(n<=p[m])return 1;
    return des(Solve(n,m-1),(ll)p[m]*Solve(n/p[m],m-1)%mod);
}
int main()
{
    init();
    int T,Case=1;
    scanf("%d",&T);
    while(T--)
    {
        ll L,R,K;
        scanf("%I64d%I64d%I64d",&L,&R,&K);
        printf("Case #%d: ",Case++);
        if(!check(K))printf("0\n");
        else if(K>M)
        {
            if(K>=L&&K<=R)printf("%d\n",K%mod);
            else printf("0\n");
        }
        else
        {
            int pos=1;
            while(p[pos]<K)pos++;
            pos--;
            int ans=des(K*Solve(R/K,pos)%mod,K*Solve((L-1)/K,pos)%mod);
            printf("%d\n",ans);
        }
    }
    return 0;
}