ZOJ 3547 & HDU 4059

/****************************************************************************************************
    大连赛区现场赛的银牌题吧(罚时少的话三题银牌)，尼玛，搞了半天没搞出来，最后发现原来是那个公式算的时候会溢出。。。
解法就是分解n之后得到n的不同素因子的个数，比如n = 12, 12 -> 2, 3，那么显然先把1^4+2^4+...+12^4加上，然后再减掉
2^4+4^4+6^4+...+12^4和3^4+6^4+...+12^4，最后再把6^4+12^4加回去，因为6和12被加了1次减了2次，这个就是dfs(后来改
位运算模拟了)+很水的容斥，不过阴险的就是那个公式比较难推，还容易溢出。。。
****************************************************************************************************/
#include <iostream>
#include <functional>
#include <algorithm>
//#include <hash_map>		//Visual C++ lib
#include <complex>
#include <cstdlib>
#include <cstring>
#include <fstream>
#include <sstream>
#include <utility>
#include <bitset>
#include <cctype>
#include <cstdio>
#include <limits>
#include <memory>
#include <string>
#include <vector>
#include <cmath>
#include <ctime>
#include <queue>
#include <stack>
#include <list>
#include <map>
#include <set>
using namespace std;

#define LOWBIT(x) ( (x) & ( (x) ^ ( (x) - 1 ) ) )
#define CLR(x, k) memset((x), (k), sizeof(x))
#define CPY(t, s) memcpy((t), (s), sizeof(s))
#define LEN(s) static_cast<int>( strlen((s)) )

typedef long long LL;		//GNU C++
typedef unsigned long long ULL;
//typedef __int64 LL;		//Visual C++ 2010
//typedef unsigned __int64 ULL;
typedef double LF;

typedef pair<int, int> PII;
typedef pair<LL, LL> PLL;
typedef pair<double, double> PDD;

//typedef hash_map<int, int>::iterator HMI;
typedef map<int, int>::iterator MI;
typedef vector<int>::iterator VI;
typedef list<int>::iterator LI;
typedef set<int>::iterator SI;

const int INF_INT = 0x3f3f3f3f;
const LL INF_LL = 0x7fffffffffffffffLL;		//15f
const double oo = 10e9;
const double eps = 10e-9;
const double PI = acos(-1.0);

const int MAXN = 10004;
const int MAXM = 10004;
const LL MOD = 1000000007;

int test;
LL n, rev;
bool is[MAXN];
vector<LL>prm, vt;

LL ext_gcd(LL a, LL b, LL &x, LL &y)
{
	if (!b)
	{
		x = 1;
		y = 0;
		return a;
	}
	LL d = ext_gcd(b, a % b, x, y);
	LL t = x;
	x = y;
	y = t - a / b * y;
	return d;
}
LL pow_mod(LL a, LL k)
{
	LL r = 1;
	a %= MOD;
	while (k)
	{
		if (k & 0x1)
		{
			r = (r * a) % MOD;
		}
		a = (a * a) % MOD;
		k >>= 1;
	}
	return r;
}
void preprocess()
{
	memset(is, true, sizeof(is));
	is[0] = is[1] = false;
	prm.clear();
	prm.push_back(2);
	for (int i = 4; i < MAXN; i += 2)
	{
		is[i] = false;
	}
	int ind;
	for (ind = 3; ind * ind < MAXN; ind += 2)
	{
		if (is[ind])
		{
			prm.push_back(ind);
			for (int j = ind * ind, k = ind * 2; j < MAXN; j += k)
			{
				is[j] = false;
			}
		}
	}
	while (ind < MAXN)
	{
		if (is[ind])
		{
			prm.push_back(ind);
		}
		ind += 2;
	}
	LL ty;
	ext_gcd(30, MOD, rev, ty);
	if (rev < 0)
	{
		rev += MOD;
	}
	return ;
}
LL fun(LL x)
{
	LL a = x % MOD;
	LL b = (x + 1) % MOD;
	LL c = ( (2 * x) % MOD + 1) % MOD;
	LL d = ( ( ( 3 * pow_mod(x, 2) ) % MOD ) + ( (3 * x) % MOD ) - 1 ) % MOD;
	LL e = (a * b) % MOD;
	LL f = (e * c) % MOD;
	LL g = (f * d) % MOD;
	return (g * rev) % MOD;
}
void ace()
{
	LL buf, res;
	int cnt;
	for (cin >> test; test--; )
	{
		cin >> n;
		if (1 == n)
		{
			puts("0");
			continue ;
		}
		vt.clear();
		buf = n;
		for (int i = 0; i < static_cast<int>( prm.size() ) && buf != 1; ++i)
		{
			if (0 == buf % prm[i])
			{
				vt.push_back( prm[i] );
				while (0 == buf % prm[i])
				{
					buf /= prm[i];
				}
			}
		}
		if (buf != 1)
		{
			vt.push_back(buf);
		}
		cnt = static_cast<int>( vt.size() );
		res = fun(n);
		for (int ind = 1; ind < (1 << cnt); ++ind)
		{
			int token = ind;
			LL tmp = 1;
			int t = 0;
			for (int j = 0; j < cnt; ++j, token >>= 1)
			{
				if (token & 0x1)
				{
					tmp *= vt[j];
					++t;
				}
			}
			if (t & 0x1)
			{
				res = ( res + ( MOD - ( pow_mod(tmp, 4) * fun(n / tmp) % MOD ) ) ) % MOD;
			}
			else
			{
				res = ( res + pow_mod(tmp, 4) * fun(n / tmp) ) % MOD;
			}
		}
		cout << res << endl;
	}
	return ;
}
int main()
{
	preprocess();
	ace();
	return 0;
}