http://acm.hdu.edu.cn/showproblem.php?
pid=4059
现场赛中通过率挺高的一道题 可是容斥原理不怎么会。。
參考了http://blog.csdn.net/acm_cxlove/article/details/7434864
1、求逆元 p=1e9+7是素数。所以由 a^(p-1)%p同余于1 可得a%p的逆元为a^(p-2)
2、segma(i^k)都能够通过推导得到求和公式 详见http://blog.csdn.net/acm_cxlove/article/details/7434864
3、容斥原理。还在恶补中 代码写的挺美丽
#include <cstdio>
#include <cstring>
#include <iostream>
#include <algorithm>
#include <set>
#include <cmath>
#include <vector>
using namespace std;
#define ll long long
#define IN(s) freopen(s,"r",stdin)
#define OUT(s) freopen(s,"w",stdout) const int MOD = 1000000007;
const int N = 10005;
const int M = 10050;
ll n,thr;//thr 30的逆元
vector<int>fac; bool is[N]; int prm[M];
int getprm(int n){
int i, j, k = 0;
int s, e = (int)(sqrt(0.0 + n) + 1);
memset(is, 1, sizeof(is));
prm[k++] = 2; is[0] = is[1] = 0;
for(i = 4; i < n; i += 2) is[i] = 0;
for(i = 3; i < e; i += 2) if(is[i]) {
prm[k++] = i;
for(s = i * 2, j = i * i; j < n; j += s)
is[j] = 0;
// 由于j是奇数,所以+奇数i后是偶数,不必处理!
}
for( ; i < n; i += 2) if(is[i]) prm[k++] = i;
return k; // 返回素数的个数
} ll qmod(ll x,ll t)
{
ll ret=1,base=x;
while(t)
{
if(t&1)ret=(ret*base)%MOD;
base=(base*base)%MOD;
t/=2;
}
return ret;
} ll sum(ll x)
{
ll ret=1;
ret=(ret*x)%MOD;
ret=(ret*(x+1))%MOD;
ret=(ret*((2*x+1)%MOD))%MOD;
ret=(((3*x*x)%MOD+(3*x)%MOD-1+MOD)%MOD*ret)%MOD;
return (ret*thr)%MOD;
} inline ll four(ll x)
{
return (((x%MOD)*x%MOD)*x%MOD)*x%MOD;
} ll dfs(int cur, ll tmp)//容斥原理
{
ll ret=0,f;
for(int i=cur;i<fac.size();i++)
{
f=fac[i];
ret=(ret+(sum(tmp/f)*four(f))%MOD)%MOD;
ret=( (ret-dfs(i+1,tmp/f)*four(f))%MOD+MOD )%MOD;
}
return ret%MOD;
} int main()
{
//IN("hdu4059.txt");
int ncase;
ll s1,s2;
scanf("%d",&ncase);
int prmnum=getprm(N-1);
thr=qmod(30,MOD-2);
while(ncase--)
{
scanf("%I64d",&n);
//int sn=(int)sqrt(n);
fac.clear();
ll tmp=n;
for(int i=0;i<prmnum && prm[i]<=tmp;i++)
{
if(tmp%prm[i] == 0)
{
fac.push_back(prm[i]);
while(tmp%prm[i] == 0)
tmp/=prm[i];
}
}
//while(tmp%prm[i] == 0)fac.push_back(prm[i]),tmp/=prm[i];
if(tmp!=1)fac.push_back(tmp);
//cout << "FUck= " << sum(n) << endl;
printf("%I64d\n",( (sum(n)- dfs(0,n)+MOD)%MOD + MOD)%MOD );
}
return 0;
}